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Trends in the Development of Liquid 
Propellant Rocket Engines 


E. G. D. ANDREWS, A.M.I.Mech.E., and A. W. T. MOTTRAM, B.Sc. 
(Bristol Siddeley Engines Limited) 


SUMMARY: The authors predict that rocket engines will fall into two categories. 
In both, combustion pressures will approach 200 atmospheres to attain optimum 
vehicle performance, balancing increased pumping losses and engine mass 
against increased expansion efficiency and reduced dissociation. For general 
purposes, engines will be intrinsically convenient: storable propellants will 
give specific impulses around 270 Ib. sec./lb., and engine specific weights will 
be between 0-02 and 0-03. For exacting applications, the greatest overall 
convenience may necessitate using unpleasant propellants so as to obtain the 
highest performance. Immediate research into the use and handling of these 
is essential if future engine development is not to be hampered. 


| 1. Introduction 


The design of liquid propellant rocket engines is not yet a complete science, 
and many aspects of their design can still be resolved only by extensive develop- 
ment of actual components. Nevertheless, a clear pattern is beginning to emerge, 
and recent engines manufactured by widely scattered organisations have many 
) features in common. An opinion that the art is reaching maturity might well be 
supported by the fact that many of these lend themselves fairly readily to the 
solution of problems posed by more exacting requirements. This seems to be 
an appropriate time, therefore, to determine the objectives to which our efforts 
should be directed in the future. 


The relative importance of various engine characteristics in terms of overall 
vehicle performance has already been discussed by a number of authorities: *. 
The authors intend in this paper to consider the engine design rather than the 
vehicle performance; in this context the more important considerations are the 
attainment of 

(a) the highest ‘specific impulse, 

(b) the highest propellant density, 

(c) the lowest wet engine mass, 


(d) the maximum operational convenience. 
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In any particular application there is an equivalence between specific impulse, 
propellant density and engine mass, this equivalence being heavily weighted jn 
favour of high specific impulse. The wet engine mass has the least influence, but 
it is necessary to remember that when the limit has been reached in specific 
impulse and propellant density it is the only factor amenable to improvement. 


There is no such equivalence in terms of operational convenience, which in 
the opinion of the authors may be the most important consideration. 


2. The Highest Specific Impulse 


The true dimension of specific impulse is velocity, and if consistent units 
were used the value obtained for this parameter would also be numerically 
identical with the mean velocity of the gases exhausted from the engine. In 
theory, the highest gas velocity is attained when the heat of reaction of the 
propellants is a maximum and is totally converted into kinetic energy, that is 


where H is the heat of reaction and J is the mechanical equivalent of heat. 
Calculated values for some propellant combinations are shown in Table I. 


These figures represent the theoretical limit of chemical rocket performance 
and are not attainable in practice, partly because of dissociation and_ partly 
because complete expansion of the exhaust gases cannot be achieved usefully. 


At the pressures commonly used in rocket engines at the present time there 
is no great difference (about ten per cent) between the performance of the two 
bipropellant systems liquid oxygen/kerosine and hydrogen peroxide/kerosine. In 
considering future developments, however, it is important to note that approx- 
imately one third of the chemical energy of the liquid oxygen/kerosine combina- 
tion is lost through dissociation whereas, owing to the large quantity of water 


TABLE I 
PERFORMANCE LIMIT FOR VARIOUS PROPELLANT COMBINATIONS 


: Mixture S.I.max 
Oxidant Fuel Ratio (Ib. sec. /1b.) 
Oxygen Hydrogen 8:1 574°8 
Oxygen Methane 4:1 481 
Oxygen Kerosine 3°43:1 466°5 
Oxygen Hydrazine 449 
Fluorine Hydrogen 9:1 900 
Fluorine Hydrazine | 707 
$5% H.T-.P. Hydrogen 20:1 434-9 
$5°% Kerosine 376 
100% H.T.P. Hydrogen 1 | 470 
100% H.T.P. Kerosine T2821 397°5 
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LIQUID PROPELLANT ROCKET ENGINES 


present in the decomposition products of hydrogen peroxide, the combustion 
temperature with kerosine is so low as to make dissociation effects less significant. 
Cf the remaining energy from either combination, approximately half is lost through 
incomplete expansion through the nozzle. 


To suppress dissociation it is essential to work with very high combustion 
pressures. To obtain high expansion efficiency it is essential to use very high 
expansion ratios and for operation near the ground at any rate this requires very 
high combustion pressures also. However the propellants have to be pumped into 
the combustion chamber (pressurised tanks and solid charges are not considered 
here) and high combustion pressures presuppose high pumping horsepowers. 
When considering optimum performance, therefore, it is necessary to determine 
the cycle efficiency, a practice which is almost unknown in the literature—Malina® 
making the only reference known to the authors. 


In a regenerative cycle“’, the highest efficiency is attained with the highest 
possible combustion pressure, but in a cycle in which the turbine exhausts to 
atmosphere, the pumping losses are higher and there is an optimum combustion 
pressure®?, 


Rough calculations made in the absence of precise data indicate that optimum 
combustion pressures are nearly an order higher than present practice, say 3000 
lb./in.? or so, depending on the propellants. 


The difficulties of working with combustion pressures of this order are con- 
siderable. In the chamber itself the temperatures will be higher because of the 
effect of pressure on dissociation, the density will be higher and the velocity will 
be higher. Although the volume of the combustion chamber can probably be 
reduced roughly in inverse proportion to the square of the pressure, and hence 
the surface area to be cooled is reduced while the propellant flow available to 
cool it is increased per unit area, the increase in local heat transfer rate necessary 
to keep the wall cool is nevertheless large. Present cooling techniques may go 
some way towards solving these problems, but at the moment this seems doubtful, 
and completly new techniques may have to be developed. 


Reference has been made to the reduction in combustion chamber volume: 
this probably means that the mass flow of propellants through it per unit cross- 
sectional area will be an order above the present levels, which may introduce new 
combustion problems. 


At present the efficiencies of rocket components other than the combustion 
chamber are of small importance because the propellants expended in driving 
them are a negligible percentage of the total mass flow. Severe limitations on 
their efficiency are accepted for the sake of lightness and simplicity. Now pump 
inefficiency results in an increase in propellant temperature and is not in itself serious, 
but because the propellants used to drive the turbine are also used inefficiently 
the efficiencies of the pumps and turbine are important in an optimum cycle. 
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The problem of using as much as possible of the energy of the propellants 
that drive the turbine is similar to the problem of obtaining the highest gas 
velocity from the combustion chamber, and one may assume that the gas generator 
pressure will be similar to the combustion chamber pressure. Because it is 
normally necessary to obtain maximum turbine power at sea level, however, there 
is no point in attempting to use the same pressure ratio as in the combustion 
chamber nozzle which is usually over-expanded at sea level. Indeed to make 
most use of the residual energy in the exhaust gases it may be desirable to main- 
tain at least a critical pressure ratio across the turbine exhaust, and to fit an 
expansion nozzle to obtain benefit at great altitudes. Even so the overall pressure 
ratio across the turbine is very great and cannot be employed in a single stage 
with reasonable efficiency, and at least two compounded stages must be used. 


There are four principal reasons why a regenerative cycle does not appear 
to be suitable for this application. The first is that it would be prohibitively 
heavy to enclose the turbine discs in a casing capable of withstanding the combus- 
tion chamber pressure (with a two-stage independent turbine, the casing pressure 
of the high pressure stage is only a few per cent of the combustion chamber 
pressure because of the pressure drop across the inlet nozzles). The second is that 
the problem of sealing the turbine shaft against a pressure of several hundred 
atmospheres appears to be impracticable whereas with the system proposed the 
rotor assembly can be overhung from the low pressure stage and the pressure to 
be sealed need be no more than about 30 Ib./in.2, The third reason is that a re- 
generative system imposes restrictions on the design of the combustion chamber 
and may not provide optimum conditions for obtaining high combustion efficiency. 
In practice it may be convenient to separate the turbo-pump unit from the combus- 
tion chamber in order to reduce the inertia of the parts that must be swivelled to 
control the vehicle: this is the fourth reason. 


With conventional combustion chamber pressures, the two objections to the 
use of very high nozzle expansion ratios are the size of the nozzle required and 
the fact that the improved performance is not obtained at low altitudes. Using 
very high combustion chamber pressures, however, the cross-sectional area of the 
throat is reduced, roughly proportionally to the pressure; also, because of the in- 
creased gas velocity, the area at the end of the expansion cone is slightly reduced for 
the same back pressure. A pressure ratio ten times as high can thus be obtained 
with a nozzle very little longer than conventional nozzles, and an improvement in 
specific impulse of some 20 per cent may be obtainable from this cause alone. 
For better performance at high altitudes a greater expansion ratio is desirable, 
and because there are grounds for thinking that over-expansion in the nozzle 
at low altitudes can be suppressed by suitable design, it is probable that the full 
benefit of large nozzles can be realised even for ground launched vehicles. Using 
high combustion chamber pressures, the improvement in specific impulse obtained 
with increasing altitude is less than with current combustion chamber pressures. 


The use of high energy propellants is, of course, subject to the same general 
design considerations. In the circumstances prevailing in the United Kingdom, 
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TABLE II 
COMPARISON OF KEROSINE AND ALUMINIUM AS FUELS FOR USE 
WITH HYDROGEN PEROXIDE 


Heat of Combustion 
B.t.u./Ib. fuel | B.t.u./Ib. mixture 

Kerosine 19,000 | 4,285 

Aluminium 10,300 | 5,450 


the authors do not see much prospect of high energy oxidants being used in the 
near future. Quite apart from any technical considerations, the production of 
fluorine or ozone would have to be increased by several orders over the present 
levels. By the time that this amount of capital expenditure becomes worthwhile it 
is likely that the development of nuclear power will have reached the point where 
its application to rocket propulsion becomes feasible. In this respect, a propulsion 
system intermediate between the chemical and nuclear rockets seems to offer 
attractive possibilities. The difficulties associated with the transfer of heat at 
high temperatures from a nuclear reactor to the working fluid are formidable, 
and early nuclear rockets are expected to operate with gas temperatures well 
below those theoretically desirable. A chemical reaction subsequent to preheat 
of the propellants in an atomic reactor may offer a significant improvement in 
performance, while operating with the reactor temperature sufficiently low as to 
minimise the risk of radioactive contamination of the exhaust due to structural 
failure of the pile. The nuclear reactor could also provide power to drive the 
propellant pumps. 


More proximately, the performance of chemical rockets can be increased 
by the use of higher energy fuels which may be already available in the United 
Kingdom or could be produced with less capital outlay than that required for 
high energy oxidants. Typical of these are metals such as aluminium and 
hydrazine, which is available as the hydrate in tonnage quantities. Aluminium 
metal, possibly in the form of a slurry, may be of particular benefit in hydrogen 
peroxide engines’. There are two reasons for this. First, although the heat of 
combustion of aluminium is low by comparison with that of kerosine, the reverse 
is true if the oxygen used in combustion is also taken into account (see Table II); 
thus, although aluminium has little merit as a heat source in an air-breathing 
engine, it has merit in a rocket engine. Secondly, the specific heat of aluminium 
oxide is low so that a large proportion of the heat of combustion is still available 
at the reaction temperature of hydrogen peroxide/kerosine. Because the large 
amount of water present has a relatively low molecular weight, the resuiting in- 
crease in temperature may more than compensate for the inability of the 
aluminium oxide to do expansive work. 


The choice is not limited to aluminium. There are numerous metals and 
alloys available, some of which have low melting points and high boiling points. 
They are excellent as coolants. 
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Ficure 1. Performance index /d" for various propellant combinations. 
Expansion from 500 lb./in.2 abs. to 14-7 Ib./in.? abs. 


3. The Highest Propellant Density 


The major importance of propellant density lies in its effect on the design of 
the complete rocket-propelled vehicle, where a large proportion of the “all-burnt” 
mass is taken up in the structure of the tanks and where aerodynamic drag may 
assume significant proportions. The relation between density and specific impulse is 
substantially fixed by the physical and chemical properties of the propellants and 
the question is one of choice. For this purpose the particular specification of the 
vehicle must be considered; there is no hard and fast rule applicable to all cases. 
Rough calculations show, however, that the specific impulse is between five and 
ten times as important as the density. A performance index Jd", where / is the 
specific impulse and d the relative density, is plotted in Fig. 1. 7 is an arbitrary 
exponent which signifies the relative importance of d. The mixture ratio is 
assumed constant and the curves do not therefore represent optimum conditions. 
For this reason, the performance of hydrogen/oxidant combinations is not shown: 
because of the low density of liquid hydrogen, mixture ratio is a critical variable. 


The effect of propellant density on the engine design is of more concern in 
the present context. In the light of the high chamber pressures which will be 
used, the power demanded by the propellant pumps will be of prime importance. 
For a given mass flow and supply pressure, the useful work done in pumping the 
propellants will be inversely proportional to the density. An increase in pro- 
pellant density will therefore effect a corresponding decrease in turbo-pump 
power requirements. By way of illustration, the mean propellant density of 
hydrogen peroxide/kerosine is some 29 per cent greater than that of liquid 
oxygen/kerosine. This will have a significant ‘effect on the size and weight of the 


204 The Aeronautical Quarterly 


| 
t 
d 
a 
as 
ti 
4 
t 
t 
a 
Vv 
i 
4 


erly 


LIQUID PROPELLANT ROCKET ENGINES 
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FicurE 2. Relative importance of mass and specific impulse. 


turbine. The pumps themselves will be reduced in size because of the higher 
density, and also to some extent because it should be possible to run them at 
higher rotational speeds, which will benefit the turbine also. Generally speaking, 
the propellant feed system will be the more compact the higher the density. This 
applies to pipes, valves and coolant passages as well as to the turbo-pump 
assembly. 


The weight reduction that may be expected to accrue from these considera- 
tions will be offset to some extent by the fact that the optimum combustion 
pressure will be somewhat higher for the denser propellants. For this reason, it 
is difficult to determine the overall effect on engine weight of variations in density, 
without a design study beyond the scope of this paper. 


4. The Lowest Wet Engine Mass 


In rocket engineering the term weight is used loosely and it is nearly always 
mass that is implied. One exception to this is during a static launch when the 
thrust must exceed the weight of the loaded vehicle. We will therefore consider 
the means by which the mass of the engine may be made as small as possible 
at the time when the flame is extinguished. This distinction between wet engine 
mass and dry engine mass is important because the wet engine mass determines 
vehicle performance. 


Consequent upon an increase in combustion pressure to increase the specific 
impulse, there will probably be some mass penalty. The effect of this on vehicle 
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performance is to detract from the benefits resulting from increased specific 
impulse, and in some conditions there may be no overall improvement. The 
relative importance of mass and specific impulse depends on the factors that are 
considered constant, but not all the combinations of these factors are capable of 
generalised presentation. If thrust and total impulse are considered to be 
constant, however, which is a realistic case in that it does not require complete 
redesign of the vehicle, the change of specific impulse which just offsets a change 
of “all-burnt” vehicle mass can be shown for various single-stage vehicles as in 
Fig. 2. In this figure the “original mass ratio” is that obtaining before the 
proposed changes are incorporated. The wet engine mass in current practice is 
usually from one to two per cent of the gross mass at launch. Taking a vehicle 
with a mass ratio of ten, an increase in specific impulse of ten per cent permits 
the wet engine mass to be approximately doubled before the vehicle performance 
is impaired. 


It is for this reason that raising the combustion pressure to improve the 
specific impulse appears to be justified. Moreover the increase in engine mass 
resulting from these increased pressures is probably less than one might suppose 
at first sight. 


Considering first the combustion chamber, the volume of the combustion zone 
is proportional to the time required for mixing and reaction. Because the 
combustion temperature and the concentration of the reactants are both increased 
as the chamber pressure is increased, the required time of stay of the gases in the 
combustion zone is reduced. Now for an ideal pressure vessel, its mass is pro- 
portional to the product of its volume and the design pressure: if, as happens in 
this case, the volume is reduced more rapidly than the pressure is increased, the 
mass of the combustion chamber itself is actually reduced. 


With higher chamber pressures, increased expansion down to the same 
absolute back pressure is achieved by a short length of additional nozzle just 
downstream of the new throat; this is because the throat area is almost inversely 
proportional to the combustion pressure. Because the gas velocity is increased, 
the exhaust end of the expansion nozzle is slightly smaller in diameter at the same 
absolute pressure when the chamber pressure is raised. 


One may safely say, therefore, that increase in chamber pressure should not 
appreciably change the mass of the chamber. 


It is highly probable that the mass of the turbo-pump unit would be increased. 
The horsepower generated is proportional to the injection pressure of the 
propellants and it seems probable that the mass would be increased in somewhat 
similar proportion to the horsepower. Taking the mass of the turbo-pump unit in 
conventional rocket engines as one quarter of the engine mass, we may therefore 
assume that the total effect of increasing the combustion pressure by an order 
would be to increase the engine mass by a factor of about three. 


As the pressure for optimum performance is approached, the engine mass 
increases more rapidly than the specific impulse increases, and the effect of this 
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change of mass must be taken into account in determining the optimum 
chamber pressure. 


5. The Maximum Operational Convenience 


Development of rocket engines to have the highest performance is essential, 
but it must never be forgotten that achievement of this is not an end in itself: in 
every case the engine is required to propel a vehicle to perform some mission, and 
each limitation on use that must be imposed in order to run the engine reduces 
the effectiveness of the vehicle. This is so, regardless of the mission for which the 
vehicle is designed, although the emphasis may change from one mission to 
another. 


Convenience is a qualitative matter difficult to discuss without vagueness, but 
it is so important that reference to it must be made here. 


Considerations of convenience begin with raw materials and end with com- 
pletion of the mission. Inasmuch as convenience determines the amount of 
human effort required to fulfil the mission, either in terms of safety precautions, 
servicing, training, provision of fixed installations, size of vehicle, duplication of 
vehicles or research and development, the measure of convenience is the cost of 
the programme. Too little is known at present about the influence of technical 
decisions on the total cost for this to be a precise criterion, however, and consider- 
able discretion must be exercised. 


The authors believe that no engine characteristic is more important than 
reliability. This must be qualified. Absolute reliability, in the sense that with- 
out servicing nothing whatever can go wrong with any engine at any time during 
its rated life, can be guaranteed only by an infinite development expenditure—it 
is suggested that development cost is an inverse function of the statistical unre- 
liability—even assuming that inspection techniques are adequate to ensure 
complete consistency between engines. In most cases it is cheaper to accept the 
probability of a few failures and allow for them by periodic inspection or provision 
of additional vehicles. 


There are, however, some types of failure which are acceptable and some 
which are not. Because rocket engines are normally used on one mission only, 
some deterioration may be permissible provided no serious mechanical damage 
ensues. For ground launched applications, occasional failure to start may be 
permissible. For air launched applications, such as second or third stages, 
failure to start is as undesirable as complete shutdown of the engine, but even 
this is preferable to failures which create conditions in which explosions or 
detonations can occur. Failures of the latter type are quite intolerable, and all 
engines with which the authors have been associated have had control systems 
so interlocked that failures of that type were most unlikely. Other safety features 
were incorporated also to prevent serious failure, even though the additional parts 
did lower the theoretical reliability of the engines. 
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One of the most serious causes of engine failure is corrosion, which is most 
difficult to overcome completely. Techniques have been developed to circumvent 
this problem—the Gamma engine in Black Knight owes much of its reliability to 
these—and recent development in the use of non-metallic materials appears to 
offer hopes of overcoming the problem entirely. These, with other developments, 
indicate that rocket engines will become increasingly reliable, until the user 
need adopt no check firing procedures. This will represent a marked advance in 
operational convenience. It should perhaps be noted that for many applications 
engine reliability might well be increased if the user accepted the rocket engine 
as basically a one shot device: simply because it is possible to do so, he often 
requests that the engine be capable of repeated firings without stripping, and 
more elaborate control systems are often fitted to cater for this. Long life rocket 
engines can be, and have been, developed for suitable applications®’, however. 


The oxidant used is often the cause of much inconvenience. By their nature, 
oxidants are unpleasant substances, and all oxidants have peculiar characteristics. 
The suitability of the oxidant for any particular mission may be more important 
than obtaining maximum performance. Thus for a defensive weapon the ability 
to maintain the vehicle in a state of instant readiness is the most important 
consideration, and solid propellant rockets are often used for this reason alone. 
The authors can see no reason why some liquid propellant combinations should 
not also be included in this category. Although liquid propellant rockets prob- 
ably cannot be considered as serious rivals to solid propellant rockets if the burn- 
ing time is less than six seconds, they become increasingly competitive as the 
duration is increased. They may be supplied as packaged units and serious 
development to this end is inevitable within the next few years. 


For tactical weapons, mobility may be most important, and vehicles which 
weigh but a few hundred pounds when empty, and can be transported and erected 
on the launcher pad with comparatively light equipment, may be desirable. In 
such cases the safety and stability of the propellants may be the prime requirement, 
since the vehicles must be filled and sometimes emptied under difficult conditions. 


Offensive weapons intended to act as deterrents fall into the same category 
as defensive weapons. It is difficult to justify the use of highly volatile propellants 
for any of these applications. Volatile propellants must be refrigerated contin- 
uously, and the rest of the vehicle must be heated. A serious penalty has also 
to be paid in performance for the weight of the lagging. It has been suggested 
that the time taken to fill the vehicle when required is acceptable because a 
longer period is required to prepare the guidance equipment. This may be true, 
but occasions can be foreseen in which the limitation could seriously reduce the 
effectiveness of the weapon. Oxidants which are suitable for these applications 
are oxides of nitrogen, nitric acid and hydrogen peroxide, and of these the last 
appears to have most in its favour. 


It is probable that the most convenient propellants for vehicles which launch 
satellites and space probes are those with the maximum specific impulse, and the 
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LIQUID PROPELLANT ROCKET ENGINES 


limitations in convenience imposed by those propellants, which are usually volatile 
—and, in the case of fluorine compounds, toxic—must be offset as far as possible 
by consistency in performance of the engine. 


Consistency of operation may be important also in engines for winged 
vehicles. For these, weight is of far less importance than for ballistic vehicles, 
and because the propellants constitute a much smaller proportion of the whole 
mass, the tendency to freeze the vehicle is reduced to manageable proportions. 
Volatile propellants may be acceptable, and because range is a function of 
specific impulse they may be essential. 


As far as the engine is concerned, there is not much to choose between any 
of the oxidants: reliable engines can be made to burn any of them, though some 
are more amenable than others. Thus hydrogen peroxide engines are simple but 
the choice of materials is limited, and there are some problems connected with 
cooling the chamber. Many materials may be used with liquid oxygen, but igni- 
tion requires more equipment, and cooling problems in the chamber, and freezing 
problems elsewhere, may sometimes be ineluctable. The authors have no 
experience of the development of fluorine rockets, but expect that engines could 
be made reliable for this too, once a range of suitable materials and the practical 
characteristics of fiuorine had been established by suitable research. Undoubtedly 
the extreme toxicity of fluorine would greatly hamper the rate of development of 
the engine until a technique had been found which would render it harmless 
without involving fire or explosion. The other objection to the use of fluorine, 
the toxicity of the exhaust, requires that vehicles be launched in remote places by 
remote control. These difficulties appear great, but the authors suspect that 
a certain amount of research into ways and means of overcoming them would show 
them to be not insuperable. Because of the increase in specific impulse obtainable 
by using fluorine, a satellite launcher need be only about one quarter the size of one 
burning liquid oxygen, and for such applications fluorine may perhaps be a 
convenient oxidant. 


6. Conclusions 


From this survey it appears that there are two parallel courses which must 
and will be followed, the two objectives being the best engines for general purposes 
and the best engines for very difficult, specialised applications. 


The best engines for general purposes will have sea level specific impulses 
approaching 270 Ib. sec./Ib., attained by using near-optimum cycles with storable 
propellants, and the main emphasis will be on dependability. This might be 
achieved by using hydrogen peroxide with selected fuels, perhaps containing 
metals as slurries. Considerable research and development is needed to attain 
this objective, attention being given to high pressure combustion and cooling 
problems, component efficiencies, and the stability of the propellants. It is prob- 
able that such engines would have wet weights of about two or three per cent of 
their thrust. 
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For very difficult, specialised engines, the programme may be extended 
through liquid oxygen, liquid fluorine, mixed nuclear/chemical, and finally nuclear 
engines. Probably of these, liquid fluorine and nuclear engines would be used 
mainly for second and third stages of composite vehicles, while liquid oxygen and 
mixed nuclear/chemical engines would be useful for ground launched vehicles, 
For these engines convenience in handling must probably be sacrificed to some 
extent. Again, near-optimum cycles will be used. The research and development re- 
quired will be far greater than for the conventional engines, and although 
development cannot be started until specific requirements are formulated, research 
activities cannot be delayed without making future development much more 
difficult. 
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Quasi-steady Derivatives for the Subsonic 
Flow past an Oscillating Aerofoil 
in a Porous Wind Tunnel 


D. G. DRAKE 


(Department of Applied Mathematics, University of Liverpool) 


SUMMARY: The compressible subsonic flow about an oscillating two-dimensional 
aerofoil in a wind tunnel with porous walls is considered. The porous wall is 
assumed to exhibit the property that the ratio of the normal velocity at the wall 
to the pressure drop across the wall is constant. The aerodynamic derivatives 
are found for quasi-steady pitching motion of the aerofoil. The dependence 
of the derivatives on Mach number and wall porosity is displayed graphically. 


1. Introduction 


A brief consideration of the problem of compressible subsonic flow about an 
oscillating two-dimensional aerofoil between porous walls has previously been 
given by the author. In the present paper the analysis is extended and the 
aerodynamic derivatives for pitching motion about the mid-chord of the aerofoil 
are found when the motion may be considered quasi-steady. 


The expressions for the derivatives include as special cases those of the closed 
tunnel and the open jet, corresponding to zero and infinite porosity respectively. 
The flow about an aerofoil at rest in a porous tunnel is also a special case of the 
problem and the results in this steady case are shown to be in agreement with those 
of Woods. 


The possibility of the elimination or reduction of the interference effects of 
the tunnel walls on the lift and moment forces acting on the aerofoil is discussed 
in Section 4. Graphs are presented of the derivatives against the Mach number 
for various values of the porosity of the tunnel walls. 


Consider a uniform stream of velocity U flowing past a thin aerofoil of 
chord c which lies midway between the walls of a two-dimensional porous wind 
tunnel of height b. Choose a co-ordinate system such that the x- and y-axes are 
parallel and normal to the walls respectively, with origin at the leading edge of the 
aerofoil when at zero incidence to the uniform stream. The co-ordinates x and y 
are assumed dimensionless, with the aerofoil chord c as the unit of length. Let 
the aerofoil perform simple harmonic oscillations of complex amplitude —Ua2 (x) 
with angular frequency » such that the normal velocity component is 
-U & [z(x)e™] where F denotes “real part of.” 
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The restriction | (x)|< 1 allows the linearised theory to be used and, since the 
flow about the aerofoil is irrotational, the perturbation velocity vector q satisfies 
the condition 


curlq=0. . ‘ ‘ ‘ 
Equation (1) implies the existence of a perturbation velocity potential such that 
iwt 


q=Uc grad ¢(x,y)e”. (2) 


Linearising the equations of motion and continuity gives 
and (M?—1) 9,, =0, ‘ 


respectively, where M is the free stream Mach number, k=wc/U is the frequency 
parameter, p is the pressure, p, is the free stream pressure and p, is the free stream 
density. The perturbation velocity potential # is anti-symmetrical® in y and 
thus the pressure jump across the aerofoil is given by (3) as 


{y (x) e™}, 
where y(x)= (= + ik) (x, 0+). — 


Some assumption must be made concerning the nature of the flow through the 
porous wall. As in previous papers” *»” the porous wall is assumed to exhibit 
the property that the normal velocity at the wall is proportional to the pressure 
drop across the wall. The constant of proportionality is taken as p,U/o with ¢ 
known as the “porosity.” With the further assumption that the pressure outside 
the wall is the free stream pressure, the boundary conditions at the tunnel 
walls”) are 


dy to (¢2+iko)=0, y=+b/(2c)=+y, (Say)... 


Also, the pressure along the centre line of the tunnel must be continuous except 
across the aerofoil and, since ¢ is anti-symmetrical in y, this implies 


¢,+ikg¢=0, y=0, x<0 and x>1. (8) 


The solution of (4) subject to the boundary conditions (7) and (8) is required. 
To this end the system (4), (7) and (8) may be subjected to the following Fourier 
transformation”. 
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. « 


which leads to the integral equation 


a(x)=- lim $, (x—Zz, y) y (2) dz . (10) 


y>0+ OY 
for the pressure jump function y (x), where ¢, (x, y) is the inverse transform of 


_ Acosh A (y, —|y|) + io (v +4) sinh A (y, — |y}) 
i(v+k)[Acosh Qy,)+ic V+k)sinh(ay,} (il) 


where sgn y=|y|/y, and with 


Se) - SF]. . 2) 


For quasi-steady flow, neglecting terms of O(k*), the inverse transform of (11) 
is (1) 


6, (x, y)= (y, —|yl) sen y— cikyB? 
co —a,x/B 
Co 1+ ik 
nm=0 ( y) 


anty,? E cos (a,y) — sin (ay) (4 + (13) 


where a,y, with 0 < tan-!(8/c) <= and x>0. 


The remainder of this paper is devoted to obtaining a solution of the integral 
equation (10) with , given by (13) and finding the quasi-steady derivatives for a 
particular prescribed motion of the aerofoil. 


NOTATION 
A’=2By, 
b tunnel height 
C=0-5772 . . ., Euler’s constant 
C, lift coefficient defined by (32) 
Cn moment coefficient defined by (34) 
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c chord of aerofoil 
f,,f. defined by (29) and (30) 
i imaginary unit 
k=woc/U, frequency parameter 
los laayat aerodynamic lift derivatives for pitching oscillations 
M Mach number of free stream 
Ma, Maajat aeTOdynamic moment derivatives for pitching oscillations 
pressure 
p, free stream pressure 
q perturbation velocity vector 
S(e,a) defined by (75) 
S’(e,2) defined by (98) 
t time 
U_ free stream 
x,y non-dimensional (with respect to c) Cartesian co-ordinates 
y, =b/(2c) 
z dummy variable 
a defined by (19) 
amplitude of pitching oscillations 


a 


y (x) pressure jump function defined by (6) 
¢ defined by (74) 
@ defined by (23) 
6, defined by (24) 
defined by (12) 
vy Fourier transform variable 
p density 
p, free stream density 
porosity 
o’=2cy,, porosity parameter 


@ non-dimensional (with respect to Uc) perturbation velocity 
potential 


®, Fourier transform of ¢, 
@, defined by (13) 


angular frequency 
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sgn y=|y|/y 
KR “real part of” 


a (x) first logarithmic derivative of gamma function 


da 
log I'(z) second logarithmic derivative of gamma function 


2. Solution for the Pressure Jump Across the Aerofoil 


In order to solve the integral equation (10) for y it is necessary to know the 
form of a(x), together with the expression (13) for ¢,. The form of a(x) is 
determined if the mode of oscillation of the aerofoil is prescribed. The particular 
case considered is that of an aerofoil oscillating in pitch about the mid-chord of 
the aerofoil, when 


Substituting for (x) and ¢,, from (14) and (13), into (10) leads to 


The assumption is made that the tunnel height is large compared with the chord of 
the aerofoil. Thus, terms dependent on negative powers of y, greater than two 
are neglected in (15). This assumption is usually made when calculating 
corrections to be applied to data obtained in subsonic tunnels. 


Now, in the Appendix it is shown that 
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and Euler's constant C=0°5772 ...... (20) 


and the logarithmic derivatives of the gamma function are defined by 


a 
d* 1 
and log (a)= (22) 
The functions defined by (21) and (22) are tabulated in Ref. 5. 
Introducing and 6, by z= —cos 6). (23) 
1 


and assuming that the pressure jump function y may be expressed 


y (0)=(A, + ikB,) cot + (A,+ikB,) sin & 


where the A, and B, are independent of k, the integral equation (10) becomes 


ikcos6]_ 1 f { } 
[1 ] (A, +ikB,) cot; + (Ant sin nO ¢ x 


ik (cos 8—cos 6,) 


26? 


x{§ exp [—a, (cos cos 6,)/8}x | 1+ 


The expression (25) for the pressure jump function is such that the Kutta condition 
at the trailing edge is satisfied and has the proper type of singularity at the leading 
edge. Evaluating the integrals of (26) and equating the coefficients of 
cos n0,, n> 0, leads to two infinite sets of equations for A, and B,, n=0. .. . ©. 
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These are easily solved and, to the order of approximation considered, give 
(1 
B [1 DBy, a) + i(3- «) + + pays! 
A,=- 8 G- a) - al. 


A,=9, n> y 


-[1- - + + (G- «) |}. 


B,=0, n> 


where 


oA, d? oA, oA, da 
(6? +0) da °S (B2+0%) | 4x da 


(30) 


August 1959 217 


|_| 
| 
|_| 
|_| 
|_| 
| 
|_| 
|_| 
(26) 
ition 
ding 
of 

, 
arterly 


D G. DRAKE 


Here A, and A, are given by (27). Thus, the pressure jump function y is 
6 
y + ikB,) cot +(A,+ikB,) sin 9+ ikB, sin 26 (i) 


where the A, and B, are given by (27) and (28). 


3. Aerodynamic Derivatives 
The lift coefficient C, is related to the lift acting on the aerofoil by 
1 


and, from (5), C.=4 | y (z) dz. (33) 


0 


C, (32) 


Also, the moment coefficient C,, is derived from the moment about the mid- 
chord by 


moment 
1 
and, from (5), (4-—z)y(z)dz. . ; (5) 


0 


Introducing the variable 6, which is related to z by (23), in (33) and (35) gives 
0 


0 


On substituting for y (6) from (31) into (36) and (37), and evaluating the integrals, 
the lift and moment coefficients are 


(A, +ikB,) += (A,+ikB,) . « @ 


since the terms O (k*) are neglected, it is convenient to write 


and [mat ikmaaat +O (k*)] . . (41) 
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Here |a, laojat, Ma and m,,,, are the aerodynamic derivatives and, in this quasi- 
steady case, are independent of the frequency parameter. Then, from (31), (36), 
(38) and (40), it is seen that 


Similarly, from (31), (37), (39) and (41), the moment derivatives are 


The derivatives for the closed tunnel and the open jet are obtained from 
(42){(45) on letting the porosity « be zero and tend to infinity, respectively. 


Considering the closed tunnel case (c =0), then, 


and C+ log (@) | 1 = —2]og 2. i . (47) 


Thus, the aerodynamic derivatives for the closed tunnel (J.). and 
are 
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(ad e= | -(1+ 5 )tog{ + 5(36°-1) (1+ 


] 


( 
| 
These results are in agreement with those of Jones“, 
For the open jet (¢ —> 0%) 
a= tan (8) |... =0 & 
and, using the relations 
d 
d? 


the aerodynamic derivatives for the open jet (1Ma)o and 
(Maa,at) o> are 


Qn 


léy, 28 \4py,/ 28 * 192B%y,2 96B%y,2 

192 los + + apy, + 2S * 


« {[1- ‘By, * tat 
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and (Mtaa/at) o= 26? { 1 4By, 968?y,? {[1- iby, * 


* 28 \apy,/ 28 * 1928.2 96B%y,2) ~ 


+ ~ 46? By, + 966%y,21 * 


1 M?x? M?x3 


In Figs. 1 to 4, the derivatives Jo, ma, Iioae aNd Mioae are plotted against the 
Mach number M for a tunnel height to chord ratio of five, i.e. y,=2-5, and for 
several values of the porosity parameter o’=2cy,. These include the closed tunnel 
(c’=0) and the open jet (o’=00) curves. 


The aerodynamic derivatives and (maar). for an 
unbounded stream are obtained on letting y, tend to infinity in (42)-(45), i.e. 


It is clear from Figs. 3 and 4 and equations (61) and (62) that the values of 
the derivatives Ijo,4: aNd Ma/a for the closed tunnel and the open jet differ greatly 
from the unbounded stream values. Furthermore, the large differences in magni- 
tude are not made smaller by using a porous wind tunnel. Thus, the effect of the 
tunnel walls on the derivatives Ijo,a: aNd Moa: Temains important at low frequency 
parameters. 


4. Steady Case 


The lift and moment coefficients of an aerofoil at rest in a subsonic tunnel 
with porous walls and at an incidence 2, to the undisturbed stream are given by 
(40) and (41) with the frequency parameter k equal to zero. Then, letting (C), 
and (C,,), be the steady lift and moment coefficients, 


(C) ola 
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and (Cy) 


E tan-? (5) + ( tan (3)) + (64) 


These are in agreement with the steady results of Woods”. 


An informative method of displaying the dependence of the steady lift 
coefficient (C,;), on the porosity, Mach number and the tunnel height to chord ratio 
is to plot B(C,),/(2™2,) against A’=2fy, for different values of the porosity 
parameter o’=2cy,. This is done in Fig. 5 for several values of o’. The difference 
in values of 8 (C,),/(22,) for a point on such a curve and a corresponding point 
on the curve £(C,),/(2™z,)=1 is a measure of the interference due to the tunnel 
walls on the lift coefficient. Thus, to measure interference-free steady lift 
coefficients in a porous tunnel, for a particular value of the tunnel height to chord 
ratio and over a range of Mach numbers, it is necessary for the porosity to be such 
that the curve coincides with 8 (C,),/(2™z,)=1. From (63), it is seen that this is 
the case if 


[tan-*(5) —4By, tan-*(5) + 5 =) (65) 


To the approximation considered, (65) is satisfied if 


Thus, the curve labelled o’=0-4 in Fig. 5 is seen to be nearly coincident with 
B(C,),/(2e2,)=1 and, for a porous tunnel with wall porosity given by (66), no 
correction need be applied to measured steady lift coefficients. 


Similarly, the interference on the mid-chord moment coefficient is zero if 
o’=?/48. Thus, while the lift and moment coefficients are interference-free at 
o’=%?/24 and o’=x?/48 respectively, the pressure distribution over the aerofoil 
at these porosities is changed. Further, the curve in Fig. 1 for o’ =0-4 (= 27/24) 
and the curve in Fig. 2 for o’=0-:2 (~ =?/48) lie surprisingly close to - closed 
tunnel (o’ =0) curves. 


Some insight into the magnitude of the error due to neglecting terms O (1/y,°) 
may be obtained in the steady case by including terms O(1/y,°) in the solution. 
Then, it is found that 


(C).= {1 apy, (G) + [tan + 
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The effect of including terms O(1/y,°) is shown in Fig. 5. It is seen that the 
higher order terms are probably only important for very large values of o” and 
then only for Mach numbers approaching the transonic régime. 


5. Conclusion 

The quasi-steady derivatives for the subsonic flow past an oscillating aerofoil 
in a porous wind tunnel have been found. Comparison with the unbounded stream 
results show that for low frequency oscillations the resultant effect of the inter- 
ference of the walls on the derivatives /, and m. is zero for all Mach numbers for 
particular values of the porosity, whereas the effect of the walls on the derivatives 
liajat ANG Maoat iS large for all values of the Mach number and porosity. It appears 
unlikely, at least for low frequency oscillations, that wind-tunnel wall interference 
is eliminated or substantially reduced when using porous wind tunnels. 
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Appendix 
The sums of the infinite series 
exp m=0, 1, 2, . (68) 
n=0 


where a,y,=tan~! (8/o0)+nz, are found in this Appendix. 
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(a) When m=0 the sum of the infinite series 


i 48 


{-4 —* a}. (69) 


is required, neglecting terms O(1/y,*). This is a geometric series and it is seen that 


Sew {-4 (x— =exp { - 2By, J 2sinh [((x—z) =/(2By,)] 


and, neglecting terms O (1/y,”), 


Me 


(71) 
(5) When m=1, the infinite series 
ae exp = B a} a," 


is to be summed. 


Tv o 


(x—z)x 
and s= —_—_,, 4 ‘ (74) 
By, 
the sum of the series — exp {—<¢(a+n)}= (<, 2) (75) 
(a+n) 
1 

is required, where a= we exp {—< (a+n) }. (76) 


Now, differentiating (76) with respect to ¢, 


co 
exp {—<«(z+n) } 
n=0 
= —exp (— ae) 
] 1 ] 


neglecting terms O 
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On integrating (77), 
. where C (a) is a function independent of «. 
Now, differentiating (76) and (78) with respect to a gives 
os 1 
) S (e, a) exp { (a+n)} (79) 
Equating the expressions (79) and (80) for 0S/da and letting ¢ tend to zero gives 
dC oo 1 
da n=0 
and, in the usual notation, 
) dC da? 
On integrating (82), 
d 
3) C(a)=-— log I‘ (a)+K,, 4 ‘ (83) 
where K, is independent of z. Now, by (77), 
oS 1) 
—exp (—<) 1—exp(—s) > (84) 
) and integrating, 
1)=—log {1-exp(—s)}+K, . (85) 
where K, is aconstant. Also, by (76), 
e>0 
and letting <=—>0o0 in (85) and using (86) leads to 
Thus, using (87), the expression (85) for S (<, 1) becomes 
) S (c, 1)=—log {1—exp (—<)} 
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Comparing (88) with the result obtained on putting 2=1 in (78) gives 


C (1)=0. (89) 
Now, the first logarithmic derivative of the gamma function is given by“ 
d 1 ©& a 
et n=, n(a+n) (90) 
where C is Euler’s constant and, putting z=1, 
~ tog 91 
dz 4 = ( ) 
Thus, by (83), (89) and (91), 
K,=-C (92) 
d 
and C(z)=-— log (2)—C (93) 
i. & ( 1 1 ) 
=— 
n+a en (94) 
Substituting for C(a) into (78) gives 
a?—@ 
and the sum of the infinite series (72) becomes 
exp {-4 *(x- a} Gn log By, log I’ (a) -C+ 
a? — } 
a— x- z)~ (96) 
where terms O (1/y,) are neglected. 
(c) The sum of the infinite series 
co 
{ - a} (97) 
may easily be obtained from the expression (95) for S («, <). 
Letting S’ a= exp {—¢(a+n)} - G (98) 
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and differentiating (98) with respect to «, 


Then, on integrating (99), 


S’ loge—e+ log I’ (a) e+Ce— e? 34K,, 


(100) 
where K, is independent of <:. Letting <—>0 in (100) and using the relation 
liimeloge=0 . 2 (101) 
«>0 
K 1 
gives (a+n)* 
Thus, by (100) and (102), 
d @ I d? 
S’ (, a)=e log e—et+ qq os + daz 08 
(103) 


and the sum of the infinite series (97) becomes 


By, By, By, 


(x— 


By, log T'(a)+C | (e-5)+ lero}, (104) 


d 


neglecting terms O (1/y,). 
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Plasto-Elastic Stress Distributions in Lugs 


H. FESSLER, M.Sc.(Eng.), Ph.D., A.M.I.Mech.E. and D. J. HAINES, B.Sc. 
(The University of Nottingham) 


SumMaryY: A photoelastic layer was bonded to the face of aluminium alloy lugs 
and the strain distribution in the lugs was determined from the photoelastic 
fringe patterns which were recorded in a reflection-type polariscope. The semi- 
circularly-ended lugs were loaded through neat fitting hardened steel pins. 
Assuming the maximum shear stress criterion of yield, the region of yielding 
under different loads was determined from the fringe patterns. Yield contours 
are shown for each of the four ratios of (hole diameter)/(lug width) tested. 
The greatest shear strains in the lugs were related to the applied loads. The 
extent of yield across the section normal to the direction of loading was 
measured for different loads and it was found that the progress of yield across 
this section is independent of the (hole diameter)/(lug width) ratio. The 
residual stresses on this section were also estimated. Complete yielding across 
the lug was related to the properties of the material and to the fracture strength 
of the lugs. Local yielding due to surface irregularities remained local, 
showing that “bedding in” does not weaken the component under static loading 
conditions. 


1. Introduction 


Static tests of lugs have been carried out by many investigators (see, for 
example, Refs. 1 and 2), in attempts to relate the nominal stresses and the strength 
of the material to the permanent deformation of the loading hole and to the failure 
load of the lug. The available information on fatigue strength of the lugs has been 
analysed by Heywood"). Linge'*’ used a photoelastic reflection technique to study 
the effect of pin clearances in lugs. 


This paper is a study of the strain distributions in lugs. The regions of plastic 
deformation were determined by assuming the maximum shear stress criterion of 


yield to apply. 


NOTATION 
b_ width of lug 
d diameter of hole 
t, thickness of layer 


tn thickness of metal 
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STRESS DISTRIBUTIONS IN LUGS 


y force to yield net normal section in pure tension 
force 


load to fracture lug 


~ force to fracture net normal section in pure tension 


y — lowest force to cause yield across one limb of a lug 
force to yield net section in pure tension 


Definitions 


NORMAL SECTION of a lug is the section normal to the direction of loading 
passing through the centre of the pin hole. 


YIELD FRONT is the boundary between elastic and plastic regions. 


2. Experimental Technique 


A photoelastic reflection method “ * *) was used. This consists in bonding a 
thin layer of photoelastic material to the metal surface and observing the fringe 
pattern in the layer caused by the strains applied to the component. These fringes 
are proportional to the difference of principal strains or maximum shear strain at 


any point. 


In this invesugation lugs made of aluminium alloy D.T.D. 646B (Ref. 8) were 
tested by bonding a layer of the epoxy resin Araldite Casting Resin B (Ref. 9) to 
the surface with the adhesive Araldite 102 (Ref. 10). The joints withstood strains 
of more than four per cent. Because the stiffness of the layer was never greater 
than 1/80 of the elastic stiffness of the metal, the layer did not appreciably affect 
the behaviour of the metal for moderate strains and thus underwent the same 
strains as the metal lug. 


Figure 1 shows the reflection polariscope and parts of the tensile testing 
machine used in the experiments. 


All lugs were manufactured from one } in. thick sheet of D.T.D. 646 B and 
had their axes parallel to the grain of the sheet. All lugs were ? in. wide and, 
except where otherwise stated, } in. thick. All surfaces were surface ground, the 
photoelastic layer being ground to the required thickness after bonding. The 
thickness of the layers, which was constant on each specimen, varied from 0-012 in. 
to 0:034 in. This variation was necessary to obtain convenient fringe orders for 
predetermined strain differences. The thickness of the adhesive varied between 
0-001 in. and 0-002 in. but was constant for each lug. The ends of the lugs were 
semi-circular, concentric with the neat-fitting, air-hardened, steel loading pins. The 
loading system used was self-aligning, thereby eliminating any applied bending 
moment on the specimens. The lug fixing was by means of a 7; in. diameter pin, 
positioned 14 in. from the loading pin. 
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FiGURE 1. The apparatus. 


To avoid errors due to creep in the adhesive, the testing time was kept to a 
minimum by photographing the fringe patterns. Most of the tests avoided 
unloading of the lugs, but it was found that repeated loading (up to 20 cycles) did 
not alter the strain distribution. 


3. Calibrations 


To evaluate the results of any plasto-elastic investigation, a knowledge of the 
stress-strain relationship of the material is essential. Conventional tensile tests 
were carried out on specimens loaded along and across the grain of the aluminium 
alloy sheet. The results are shown in Fig. 2. Young’s modulus for both specimens 
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\ 


STRESS 


nm 
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2 3 4 5 6 
DIRECT STRAIN (1 UNIT = 1% OR 0:1% 


FicurE 2. Stress-strain relation for D.T.D. 646B in pure tension. 
©—test in direction of grain. X—test perpendicular to grain. 


was 10-1 x 10° Ib./in.? and Poisson’s ratio was 0-30 for the elastic state and 0°50 for 
the plastic range. 


Araldite Casting Resin B, the material of the photoelastic layer, behaves 
elastically over the test range. Young’s modulus was found to be 0-48 x 10° Ib. /in.’, 
Poisson’s ratio 0-38, and the strain fringe coefficient 11,020 fringes/in. for unit 
Strain difference. 


The adhesive, Araldite 102, has a Young’s modulus which is never greater than 
20,000 Ib./in.?, and a strain fringe coefficient of 2,200 fringes/in. for unit strain 
difference. Because the thickness of the adhesive was only 0-001 in. to 0-002 in., it 
contributed little to the optical and mechanical behaviour of the layer, which varied 
between 0-012 in. and 0-034 in. for the specimens tested. 
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Part of the load applied to the lug is carried by the photoelastic layer and this 
introduces an error in the stress distributions. If it is assumed that the strains in 
the layer and metal are equal, because the widths are the same, the error 


e— load carried by layer 
1+ Pm/P, 1+ (Entm)/ (Ext) 


‘total load 


In the elastic range the error e was never greater than 1-2 per cent. At 44 per cent 
strain the error e is 5 per cent if the thickness of the layer is 1/6 of the metal 
thickness. 


The elastic strains at the edge of the hole in the plane normal to the direction 
of loading were measured for three hole sizes and the elastic stress concentration 
factors were compared with values interpolated from the results of Frocht and 
Hill" for square-ended lugs (Table I). The differences were attributed to the 
difference between round and square-ended lugs and to errors in both sets of 
experiments and in interpolation. 


TABLE I 
COMPARISON OF STRESS CONCENTRATION FACTORS 
Authors 2°78 2°76 2°42 
Frocht and Hill) 2°74 2°45 


4. Strain Distributions 


Figure 3 shows the fringe pattern in a typical lug after the application of the 
stated loads. It was thought most useful to express the loads as fractions of y, the 
force needed to cause yield in a uniform tensile specimen of the same material and 
the same net section, (b—d) tn». For this purpose a yield stress of 65,000 Ib./in.’ 
was decided on from the test data, shown in Fig. 2. The maximum shear stress 
criterion of yield is the justification for the assumption that yielding in the bi-axial 
stress system on the surface of a lug occurs at the same stress difference of 
65,000 Ib./in.2 This corresponds to a direct strain of 0-65 per cent and to a strain 
difference of 0-8 per cent. The strain differences are measured by the photoelastic 
layer and the thickness of this layer determines the fringe order which corresponds 
to yielding. 


The photographs in Fig. 3 show, typically, general yield originating at the 
inner edge of the lug at positions inclined at 60° to 90° from the axis of loading 
when viewed from the centre of the pin. The “yield front,” that is the boundary 
between the elastic and the plastic region, advances radially and circumferentially 
under increasing load, as may be seen by studying the position of the second 
order fringe. 


This advance is most conveniently shown by the yield contours in Fig. 4. It 
may be seen that a region “behind” the pin remains elastic. As expected, the 
outer portions of the “small pin lug” in Fig. 4(a) remain elastic up to its failure by 
shear. Fig. 4(c) shows unsymmetrical behaviour. Yield first occurs in the right 


236 The Aeronautical Quarterly 


h 
li 
le 
b 
0 
a 
S 
S 
u 
a 
A 


this 
S in 


cent 
etal 


tion 
tion 
and 
the 
of 


STRESS DISTRIBUTIONS IN LUGS 


MAXIMUM SHEAR STRAIN 
UN) 
+ 


fe) 0-2 0-4 0-6 0-8 
LoaD/ (FRACTURE LOAD) 


FiGurE 6. Maximum shear strains in all lugs. 
d/b 1/3 5/12 1/2 7/12 
Symbol © x 


hand limb and spreads across the limb at a load of 0-7, when most of the left hand 
limb is elastic. The yield front becomes slightly more symmetricai under greater 
loads. The last photograph for Fig. 4(b) was taken when straining was halted 
before failure under the ultimate load. This explains the 1-00 contour. 


Initial yielding always occurred at the loaded pin-lug interface, but the position 
of this yielding and the load at which it occurred varied greatly. This was 
attributed to local surface irregularities of the components in contact. Figure 5 
shows an example‘of such localised yield. If it occurred in regions of low nominal 
Shear stress, e.g. near the axis of loading, this yielding was confined to a small area 
until the general stress level caused widespread yield in the lug. This proves the 
common belief that local surface irregularities are not detrimental to the static 
Strength of components because the consequent “bedding in” causes only 
local yielding. 


The greatest shear strains were determined from the recorded fringe patterns 
and related to the loads causing them. 129 readings were obtained from the lugs 
tested. 101 of these lie within the scatter band shown in Fig. 6. The remaining 
28 values are shown individually. The points which lie above the scatter band 
correspond to observed local irregularities near the axis of loading. The only 
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FicureE 7. Static strength based on properties of the material. 


significant difference between the different shapes of lugs is the lower strain 
difference in the small pin lugs under very large loads. 


Figure 6 gives a guide to the ratio of the “repeated loading strength” to static 
strength of lugs. From the upper boundary of the scatter band the yield shear 
strain of 0-4 per cent corresponds to a quarter of the fracture load. If it is 
assumed that reverse yielding occurs when unloading after twice the yield strain 
has been attained on loading, this would occur at 0-8 per cent shear strain or half 
of the fracture load. Hence it may be suggested that lugs remain fully elastic if 
loaded to no more than a quarter of their static strength, yield unidirectionally if 
loaded to a half of their static strength and are likely to yield twice during every 
loading cycle if the variation in load exceeds half the static strength. 


The usual measurements of permanent deformation were made on one lug 
(d/b=1/3) after the application of a load of 2,900 Ib. (which subsequently caused 
failure). The loading pin was in single shear and its elastic deflection and its 
plastic deformation (0-0013d) caused a variation in the elongation of the pin hole. 
The elongation of the hole was 1:19 per cent at the back face and 1-08 per cent at 
the front face of the lug. This does not affect the other results presented in this 
paper, because these variations through the thickness of the lug occur at large 
strains only. 


The effect of thickness variation was studied with two thin lugs. These had a 
(limb width)/thickness ratio of (b —d)/(2t)=1-95 instead of the normal (b — d)/(2t) 
=11 for the narrowest limbs tested, but the results showed no significant 
difference from the lugs of standard thickness. 


5. Ultimate Loads 


Because it is easy to determine the ratio (U) of the ultimate load carried by a 
lug at fracture to the ultimate strength of the net section subject to pure tension, 
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FiGure 5. Example of local yielding due to surface irregularities. 


d/b=1/2, load =0°53y. 


FIGURE 11(a). 


FiGure 11(b). Fringe pattern after removal 
of load of 0-80y. 
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Fringe pattern under a load 
of O'80y. 
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FicureE 8. Safe static strength based on ultimate strength of lugs. 


an attempt was made to relate this ratio with some yield characteristic of the lugs 
tested. Figure 7 shows the mean values of U and the mean values of Y, (force to 
cause yield to spread across one limb of a lug)/(force to yield the net section in 
pure tension). It is evident that both graphs are of the same form. Values of Y 
could only be determined as the mean of the last load before, and the first load 
after yield has spread across a section. The greatest divergence of any point from 
the Y curve was estimated as +0-025Y. 


Figure 8 is approximately equivalent to a plot of Y/U, because the yield 
strength of the material is almost equal to the ultimate tensile strength 
(65,000 1b./in.? and 67,000 Ib./in.? respectively). It was assumed that for zero 
hole size the yield fracture ratio would be the same as for a uniform bar. Linear 
relationships are not claimed, but are convenient because the change of slope 
coincides with the transition from “shear failure” to “tensile failure” (see Fig. 4). 
It may be concluded from Fig. 8 that, for practicable designs of lugs, complete 
yield across a section is unlikely to occur at less than 0-6 of the force required to 
fracture the lug. The authors believe that it is inadvisable to rely on the strain- 
hardening characteristics of the material and therefore suggest that the force which 
causes yielding across the whole section of one limb of the lug should not be 
exceeded. Figure 7, or Fig. 8 in connection with the ultimate load of the lug, thus 
give the static strength of semi-circularly-ended lugs. 


6. Section Normal to Direction of Loading 


Figure 9 shows typical plots of fringe order (proportional to maximum shear 
strain) across the normal section of one of the lugs tested. Each graph was 
obtained from a photograph of the fringe pattern for a particular load. The curves 
for the small loads are similar for all sizes of lugs and show an increase of strain 
gradient towards the pin hole. The horizontal lines denote the fringe order which 
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FiGURE 9. Fringe plots across normal section of a lug. d/b=5/12. The numbers against 
the curves indicate the applied loads. 


corresponds to the strain difference which is assumed to cause yielding. In many 
cases the strain gradient becomes fairly constant across the section under large 
plastic strains. This represents a strain distribution analogous to the stress 
distribution caused by an eccentric tension. 


Fringe distance curves similar to those in Fig. 9 were plotted for fifteen lugs, 
and the position of the yield front was measured in each case. These positions 
were related to the load on the lug. Figure 10 is a record of all the results. The 
scatter band contains 120 points randomly distributed. This shows that yield 
progresses across the normal section as shown, independent of the (hole diameter)/ 
(lug width) ratio and is not affected by variations in the load and position of the 
onset of yielding. The difference between small hole lugs (d/b=1/3) and the 
other sizes tested (d/b=5/12, 1/2, 7/12) is that the small hole lugs failed in shear 
when yield had progressed 0:6 of the distance across the normal section. The other 
lugs yielded across the whole normal section before failure. 


7. Residual Strains and Stresses 


The residual fringe pattern after the load has been removed can be recorded 
and these fringes are, of course, proportional to the residual difference of principal 
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Ficure 10. Position of yield front on normal section of all lugs. 
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strains. The residual stresses cannot readily be determined from this because, 
in general, 


(a) the directions of the original and the residual stresses do not coincide, 


(b) the residual stresses act in the same directions as the last plastic strain 
increment, not in the direction of overall plastic strain, 


(c) plastic stress-strain relationships are not uniquely determined by any 
simple test. 


In spite of these theoretical objections the attempt was made to determine the 
residual stresses on a section normal to the direction of loading. 


It was assumed that the applied and residual strains and stresses in the normal 
section are parallel to the direction of loading and that the strain at right angles to 
the direction of loading is negligible (or constant). This assumption was necessary 
because it was not possible to record the directions of the residual strains reliably. 
The directions of the strains under elastic and considerable plastic deformation 
were determined for one lug. These results suggested that the foregoing assumption 
may be a fair approximation. 


Figure 11(a) shows the fringe distribution in a lug of d/b=4 under a load of 
0-80y. Figure 11(b) shows the fringe pattern after removal of that load. The 
negative fringe orders indicate a reversal of the sign of the strain difference. 
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The fringe or strain distributions across the normal section under load and 
after removal of the load, shown in Fig. 12, were plotted from Fig. 11. To obtain 
an approximate stress distribution from this information it is necessary to know the 
stress-strain relationship for the material in tension and compression and to assume 
that these relations also apply for the lug considered. From Fig. 2 it was assumed 
that yield occurred at a stress of 65,000 1b./in.* in tension and compression and 
that this stress is not exceeded (i.e. no work-hardening). With these assumptions 
it is possible to plot the stress distributions shown in Fig. 12. At any point, the 


fringe order under load which corresponds to plastic strain must be subtracted from 
| the residual fringe order to obtain the residual stresses. 


It may be seen that reverse yield has occurred on unloading. From the graph 
of residual stress, the ratio of the positive area under the graph to the negative area 
is 0°97. This indicates the overall error due to the measurements and assumptions 
made. A test on another lug gave a similar result. 


8. Conclusions 


The results should apply approximately to semi-circularly-ended lugs of any 
material which has a stress-strain relationship similar to that shown in Fig. 2. 


The greatest shear strains occurred at the loaded surface of the hole, usually 
at 60° to 85° from the axis of loading. This shear strain depends on the load only 
and is independent of the (hole diameter)/(lug width) ratio. If shear strain is 
considered a suitable fatigue criterion, Fig. 6 gives relevant information for the 
greatest possible value of fatigue strength. 


The force required to cause yielding to spread across the whole section of one 
limb was related to the ultimate load in Fig. 8. Because it is inadvisable to exceed 
this force, Fig. 8 gives the safe static strength. 


The normal section is not the most highly stressed region of a lug. The 
progress of yielding across this section (shown in Fig. 10) is independent of the 
(hole diameter)/(lug width) ratio. The residual strains after removal of a load 
of 0:80y showed extensive reverse yielding. 
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Notes on Thin Wing Theory at 
Low Supersonic Speeds 


Part I. The Design of Lifting Surfaces 


G. J. HANCOCK 
(Queen Mary College, University of London) 


SumMaRY: On the basis of the linearised theory, the integral relationship for 
the incidence distribution in terms of the velocity potential is established for 
wings with subsonic leading edges. Some analytical problems are analysed. A 
simple general numerical method is given for this design problem which com- 
pares favourably with exact linear theories. In Part II, to be published, a 
further numerical method is developed, for calculating the loading on any 
specified thin wing with subsonic leading edges, which again agrees favourably 
with exact linear theory. Both of these numerical techniques can be easily 
accommodated on desk calculating machines. 


1. Introduction 


Most of the effort devoted to wing theories in supersonic motion has been 
concentrated on the problem of determining the loading on a specified finite thin 
wing. This problem may be considered as solved for the majority of practical 
cases on the basis of the linearised theory, either by analytical or numerical 
techniques. However, the inverse problem of determining a wing surface which 
conforms to a specified loading distribution has attracted far less attention, although 
the relevant formulae are presented in Ref. 5. 


In this part of the paper the emphasis is on the design problem for wings with 
subsonic leading edges. The relevant formulae are deduced from the well-known 
lifting surface equation for the sake of completeness. Although these formulae can 
be applied analytically to any particular design consideration, the integrals become 
intractable for complicated plan forms and loading distributions; therefore a general 
numerical method is developed. This method can be applied with reasonable 
accuracy for the whole range of wing parameters (curved or straight, subsonic 
trailing and leading edges) and loadings. 


No mention is made in this paper of the type of loading required to optimise 
some aspect of performance; this is left to those who might wish to apply 
the method. 
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NOTATION 


X,y,Z 
Cartesian co-ordinates 

q_ velocity vector 

V free stream velocity 

M free stream Mach number 

u,v oblique Mach co-ordinates 
w downwash velocity 


u=f (v), v=f (uw), \ 
x=I(y) j 


c, centre chord of wing 


equations of wing leading edge 


S semi-span 
d length of mesh system 

K*=(cot 0—)/(cot 

K,?=Mc,/B 

b,=K,?/(1+K?’) 

b,=K,* (1—K*)/(1+K?’) 

b, =(K,?/2) (1—K’) 

Pp, constant, representing uniform loading 

Pi=D, (cot @+8)/M 

Ai, ine,» coefficients for numerical method 

@ velocity potential 


¢:,; the mean velocity potential between the lines 
{(i-1)d <u Sid, (j-1)d<v Sj} 


@ wing incidence 


a; ; mean wing incidence between the lines 
{i-l)d<uSid, 


=0...k 


f(x,y) equation of wing camber surface 
p density of free stream 
8 distance from edge singularity 
6 half the wing apex angle 
B=(M?-1)/? 
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THIN WING THEORY 


9, Design of Camber Surfaces 


On the basis of the linearised theory of supersonic flow, the perturbation 
velocity g (u,v, w), relative to a set of Cartesian axes, is superimposed on a steady 
supersonic stream (V, 0, 0) and satisfies the condition of irrotational motion. Hence 
a velocity potential @ exists such that q=grad@. Taking the origin of the axes of 
reference at the leading apex of a swept wing (the x-axis is in the direction of the 
main stream, the y-axis is in the spanwise direction and the z-axis completes a 
right-handed set), the physical problem of the calculation of the loading distribution 
over a thin wing is resolved mathematically by the determination of the discontinuity 
of ¢ across the plan form of the wing on the plane z=0 which satisfies the condition 
that the downwash is continuous and makes the flow tangential to the wing surface. 


That is, 


BEY, 


where z=¢ (x, y) is the equation of the wing camber surface and < (x, y) is the wing 
incidence distribution. The loading is equal to 2pV0¢ (x, y, + 0)/ 0x. 


The well-known fundamental formula” applying this principle is 


w 1) dé dn 


S 


where S is the fore Mach sector on the plane z=0 upstream of the point (x, y, 0), 
ie. for (x — €)? > 8? (y—7)”, where 6? = M?-1. 

When the wing has supersonic leading edges and the Mach sector S does not 
intersect the wing tips, w (€, 1) is known over the whole of S, being specified on the 
wing and zero ahead of the leading edge. Equation (1) can then be applied directly. 
If S intersects the wing tips it will also include some of the upwash field off the wing; 
however, this is related directly to the incidence of the wing by Evvard’s technique™, 
so equation (1) is again readily applicable. For wings with subsonic leading edges 
S will include the upwash regions ahead of the leading edges but so far it has not 
been possible to relate analytically these regions to conditions on the wing, although 
approximate methods have been devised®: *. 


In the design problem we wish to specify the loading over a given plan form 
and to calculate the incidence and camber surface of the wing. Since the loading 
distribution is defined in terms of ¢ (x, y, +0) we wish to invert equation (1) to 
obtain an integral relationship for w (x, y) in terms of ¢ (€,7,+0). In this paper the 
investigation is restricted to cases of wings with subsonic leading edges only, 
although the methods of inversion can be applied more generally. 


Consider a wing with subsonic leading edges, as shown in Fig. 1, symmetrical 
about the x-axis. 
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WAKE 
WING PLAN 
FORM 


FicurE 1. Notation. \ 


On transforming the Cartesian co-ordinates (x,y) to oblique Mach co-ordinates 
(u,v), on the plane z= +0, defined by 


M 


2p (x+ By), 


M 
n= BY)» v= 


as shown in Fig. 1, equation (1) becomes 


| 

0 v=0 
As already pointed out, w (u,v) is known only within the part of S on the wing 
surface and is unknown ahead of the leading edges. However, taking w (u, v) to be 
an integrable function throughout the whole of area S, equation (2) may be written 
as the two simultaneous equations 


u 


—7™M¢ (u,,v,)= | V(u,—u) 
0 
: (3) 
F (u, v,)= 


These are both Abel type integral equations, so that, if 


(u(Ode 
where 
k 


k 
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THIN WING THEORY 


Applying the inversion formula to both equations of (3) and combining them, 


remembering that 2 (u, v)= — v) 


¢ (u, v) dudv 


a(u,,v,)= ia 02, Vit, -2)] (4) 


This is therefore the basic design formula which gives the incidence distribution on 
the wing from a known distribution of ¢(u,v). Since loading only occurs on the 
wing surface, 0¢ (x, y,+0)/0x is zero off the wing plan form, and the discontinuity 
¢(x, y,+0) occurs only on the wing surface, together with the trailing vortex wake, 
as in subsonic wing theory. Furthermore, the discontinuity in the wake is deter- 
mined by the conditions at the trailing edge and is a function of y only. Specification 
of the loading on a wing implies specification of ¢ (u,v) on the wing and the wake, 
while # (u,v) is zero outside these regions. Therefore in this case equation (4) can 
be applied, directly the loading is stated, over the whole wing; the wing tip and 
subsonic trailing edge regions have no special significance in this case. 


It should be mentioned that equation (4) is not new but is essentially one of the 
fundamental formulae in supersonic wing theory, derived from first principles by 
Heaslet and Lomax, but it appears to be not very well known. The author thought 
that the present method of obtaining equation (4) from the familiar equation (1) 
would be the most satisfactory introduction for most readers. In any case, the 
application of equation (4) (to which the rest of this paper is devoted) is not 
investigated in Ref. 5. Goldsworthy has shown that equation (4) is applicable for 
the symmetrical problem of the design of the thickness distribution with a specified 
surface pressure distribution at zero lift when the leading edges are supersonic. 


3. Form of the Velocity Potential 

In the present design problem the loading will be assumed over a given wing 
and therefore the discontinuity *~ 0@/0x is known across the wing plan form. It is 
worth reviewing at this stage some of the main points to be included in any 
assumption of 0¢ (x, y, + 0)/ 0x. 


At the leading edge : — 


(i) If the incidence at the leading edge is finite, then 0¢/0x=O (5—"/*), where 
8 is the distance from the leading edge, except at a unique value of the 
incidence (equivalent to the ideal angle of incidence in two-dimensional 
subsonic theories) when 0¢/d0x=O (8"/*). 


(ii) If the loading is assumed to be O(1) at the leading edge, the incidence 
there will be logarithmically infinite. 
At the wing tips: — 


(i) If the incidence is finite at the wing tips (assuming that they exist), 
a¢/dx=O (8,'/?), where 5; is the distance from the wing tips. 
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(ii) If 6¢/0x=O (1) at the wing tips, the incidence there will be infinite, of 
the order 5,~'. 


At the trailing edge: — 
(i) If the trailing edge is supersonic, then 0¢/0x=O (1) for a finite incidence, 


(ii) If the trailing edge is subsonic, then 0¢/dx=O(6,'/*), where 8, is the 
distance from the trailing edge, thus satisfying the Kutta-Joukowski 
hypothesis. 


On assuming the distribution of 09 (x, y,+0)/dx, the velocity potential follows 
immediately from 


€=2 


f=1(y) 


where £=/(y) is the equation of the leading edge, so ensuring that ¢ (x, y,+0) is 
zero along this edge. It is seen that, when x lies in the wake, where 0¢/0€ is zero, 
¢ is a function of y only. 


A further consideration is now introduced, namely, if 0¢ (x, y,+0)/dy is dis- 
continuous anywhere within the wing area, while the loading 09 (x, y,+0)/dx 
remains finite, the downwash will be infinite along the discontinuity. This follows 
basically from continuity considerations since, if d¢/0x is finite and 0¢/ dy is infinite, 
then w=09/0z is infinite. There are two common examples of this type of 
discontinuity. The first occurs on the centre line of swept wings with uniform 
loading. The second occurs when the lines of constant loading (or isobars) 
are kinked. 


4. Wings with Uniform Loading 


Two simple examples are worked out to illustrate Section 3. The camber’ 


surfaces of a delta and of a swept wing with subsonic leading edges are calculated to 
give uniform loading over the wings. 


FIGURE 2. Co-ordinate system of delta wing. 
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THIN WING THEORY 


41. DELTA WING 


This wing is shown in Fig. 2 and is defined in terms of the Mach co-ordinates 
(u,v). The leading edges are v=K*u and u=K*v where K*=(cot —)/(cot 6+ 
The trailing edge x=c, is therefore u+v=Mc,/8=K,’, say. 


The condition of uniform loading is 


constant = 
Ox Po- 
On application of equation (5), 


$= Py (x—|y| cot 4), 


which becomes, on transformation to the (u, v) system, 


for y>0 6 
=p,(v—K*u) for 
where = 


The incidence distribution is found by substituting equation (6) in equation (4), 
that is 


u=u, = 4, 
u=0 v= K*u = 
u—K?v 
J 
Mp. u’ 


where u’= /(u,/v,). 


This example illustrates some of the points raised in Section 3. Since the 
loading is assumed to be finite at the leading edges (u’=K, 1/K), the incidence there 
tends to — 00 like log 8. Also the assumption that @9/dy is discontinuous across the 
wing centre line (u’ = 1) gives an infinite downwash on this line. Diagrams are given 
in Fig. 3 for the incidence distribution and camber surface of a 45° swept delta 
when K?=1/4 (M=1-2) at the spanwise positions of s/4, s/2 and 3s/4. The 
camber surface ((x, y) is derived from the standard linearised formula 
=a (x, y). 
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FIGURE 3. 


4.2. SWEPT WING WITH SUBSONIC TRAILING EDGES 


In this example the delta wing is extended to have a subsonic trailing edge 
(parallel to the leading edge), as shown in Fig. 4. 


The assumption of uniform loading is 


=constant = 
Ox 


so that 


Ss=cCy 
¢=p,(x—|ylcot#) on the wing 
Cy+|y| cot 6 > x > |y| cot 4 


x>c,+|y| cot 
=DPlo on the wake 


ly| Ss=c, 
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Hence, using the (u, v) co-ordinates, 


¢=p,(u-—K*v) for y>0 
on the wing 
=p,(v-—K*u) for y<0 
=P Co 


on the wake. 
=p,b, 


Four areas of influence are denoted in Fig. 4 and the incidence distributions are 
given for: — 


Area I :—this is given by equation (7) and is denoted by 2; (u,, 7,) 


Area II:— 
u=v-—b, 
=b, 


Mp, u, + K? 


+5(1+K?) cosh ( 
Area III:— 
usu, v=u 
Mp, @& | du v—K? b, 
= 
usu, v =(u—b,)/K? 
V(u,—u) V(r, —wv) 
u=K,?/2 v=u 
M 2 (1-K? 
—b,+K?*u, —K’K,? 


v,+u,—K,’ } 
—K cosh Pe. +2 cosh? ) 


Area IV:— 


(u,, v,)= (u,, v,)+ (u,, ay (u,, 
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THIN WING THEORY 


This example also illustrates some further points of Section 3 in the respect 
that the incidence tends to infinity at the wing tips and the trailing edge because of 
the assumption that the loading is finite at these edges. 


The incidence and camber surface distribution is given in Fig. 5 for the 
spanwise stations y=s/4, s/2, 35/4. 


5. A Numerical Technique 


In general the analytical approach of specifying the loading, substituting in 
equation (4) and evaluating the numerous integrals for the incidence becomes more 
tedious and difficult as the loadings and plan forms become more complex. There- 
fore the next step is the investigation of a numerical computation of equation (4) to 
cover all cases with equal facility and to any desired accuracy. There are several 
ways of attempting to do this. These depend on the extent to which the differentia- 
tion process is taken under the integral sign analytically before resorting to 
computation, with the idea of increasing the accuracy for a minimum of 
computation. The author investigated these possibilities but came finally to the 
conclusion that the direct numerical evaluation is the simplest approach, giving 
sufficient accuracy. This direct method is the only one discussed any further. 


First it is assumed that the loading is known either as an analytic function or a 
discrete function over the plan form. Therefore the velocity potential may be 
regarded as specified over the wing and wake. On forming a mesh across the wing 
plan form and trailing wake, composed by oblique Mach co-ordinates (u=0. d, 2d, 
....nd; v=0, d, 2d... . nd, where d is a small distance), it is assumed that the 
velocity potential ¢ is constant within each area of the mesh. Denoting by 4, ; the 
constant velocity potential enclosed within the lines {(i-1)d<uSid, 
(j-1)d<v < jd} it is assumed that ¢; ; is the value of the specified function at the 
centre point of the area, that is 


$i, [((i—4) d. G—4) d]. 


These points are illustrated in Fig. 6. 


vend 
FicurE 6. Notation. 
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Now consider the integral 


u=id v=jd 


du i, 
u=(i-1)d v=Gj-1)d 


assuming that u,=kd, v,=Id is a point on the wing affected by 9, ;. Thus the 
coefficient A (i, yx, », defined by equation (8), is finite when k >i and 1 >j, but zen 
for other conditions. 


On summation over the wing, extending equation (8), 


If a, is the average incidence over the area {(kK—1)d<u<kd, (I-l)d< 
v </d}, remembering that equation (4) states that 


M 


Although equations (9) and (10) can be combined, it is better, from a numerical 
point of view, to leave them separate. 


Once the table of A4Q, jx, has been established (depending only on the integers 
k—i, 1—j), operation (9) is easily applied on a desk machine, provided, of course, 
that the mesh is not too fine. 


FIGURE 7. 
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LEADING EDGE AREAS 
GIVEN BY EQUATION (11) 


us (k-2)d 


* usd 
AREAS AFT OF us2d 


LEADING EDGE 


GIVEN BY 
EQUATION (8) 
DING 
(4%) 
FicurE 8. Mesh pattern with leading edge FiGurE 9. Sub-division of mesh pattern 
modified for numerical evaluation. for (u,,, V;) near the leading edge. 


Refinements of this method need to be introduced for areas near the leading 
edges. The author has adopted the following procedures. 


(i) If the point (u,, v,) is not “near to the leading edge” (this point is clarified later) 
the value of ¢;, ; for an area intersected by the leading edge is taken as 


S, 
(11) 


where ;=9 d, (j—-4) a 


using the notation shown in Fig. 7. 


This means that, for a point (u,,v,) inboard of the leading edge, the leading edge is 
effectively modified, as shown in Fig. 8, with the velocity potential given by 
equation (11) at the leading edge areas and by equation (8) for those areas aft of the 
leading edge. 


The introduction of this sub-mesh introduces the need for a further table 
Ai, ix, » for small non-integral values of k—i, ]—j, depending on the shape of the 
leading edge. 


The final point concerns whether or not (u,, v;) belongs to case (i) or (ii). The 
author assumed that if u, or v is less than the distance 2d away from its nearest 
leading edge, then (u;, v,) belongs to case (ii), and if farther away it belongs to 
case (i). 


(ii) If the point (u,, v) is near the leading edge the mesh is subdivided locally with 
the velocity potential calculated on the same principles as before (i.¢., the value of 
¢ is ¢’, defined in equation (11), at the centre of the sub-area, and it is assumed that 
this value holds over the whole sub-area). This is shown in Fig. 9. 
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This technique has been applied to three examples, which are discussed in th 
following section. The full numerical computation of the first example is given 
in the Appendix. 


6. Numerical Examples 


Three examples have been worked out. These are: — 


(i) To check that the loading of a flat plate delta wing gives a unifom 
incidence distribution. 


(ii) The application of the numerical technique to the uniform loading ove 
a delta wing, and comparison with the analytical solution given in Fig, 3, 
This is designed to check the accuracy of the numerical method wha 
singularities of the incidence distribution are expected at the leading edge, 


(iii) The application to the uniform loading over a swept wing with subsonic 
leading and trailing edges and comparison with the analytical solution 
given in Fig. 4. This will check the accuracy of the numerical method 
when singularities are expected at the trailing edge. 


Example (i) 


It is well known that the velocity potential on the upper surface of a flat plate 
delta wing with 45° sweep is 


(x, y+0)= (?-y)”, 


E’ 
2/2 


where E’(8)= (1 —(1 — 6?) sin? 


0 


In this example 8 =0°6 is assumed, so that 


Va 2 


in oblique co-ordinates. 


The potential equation (12) is now assumed over a 45° swept delta wing and 
the incidence is calculated numerically by the method of Section 5. The mesh is 
specified in this case to give 10 points (or areas) at the centre span and 5 areas at 
the quarter span, as shown in Fig. 10. 


The numerical method for this example is presented more fully in the Appendix. 
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It is seen that, although the mesh is comparatively wide, the error for the 
incidence is only appreciable in the area at the leading edge; otherwise the agreemen, 
is good. However, on integration of the incidence distribution for the cambe 
surface, the effective error is less near the leading edge, so that the agreement with 
the flat plate is good over most of the wing surface. 


The incidence distribution and wing shape has been calculated at the centre 
section and the quarter-span section and is compared with the flat plate distribution 
in Fig. 11. 


Example (ii) 


In this case a uniform loading is assumed over the same delta wing as in the 
previous example, moving at the same speed. The mesh system of Fig. 10 is again 
used as the basis for the computation. 


Equation (6) states that, in this case, 


where p, is a constant. 


The downwash distribution and the camber line have been computed for 
y=s/2 and are compared with the exact theory, presented in Section 4, in Fig. 12. 


Again the error is greatest near the leading edge. but the characteristics of the 
camber line are brought out, even allowing for the infinite leading edge downwash. 


Example (iii) 


Example (ii) has been extended to the 45° swept wing with uniform loading, 
discussed analytically in Section 4 (leading and trailing edges parallel, aspect 
ratio 2). The velocity potential is again given by equation (13), remembering that 
it remains constant aft of the trailing edge. The comparison of theory and 
numerical approximation is shown in Fig. 13, for both the downwash distribution 
and the camber line at the section y=s/2. 


The regions of error for this example are confined to both the leading and 
trailing edges, although the overall characteristics of the wing are obtained. The 
conditions at the trailing edge do not satisfy the Kutta-Joukowski hypothesis, but 
more orthodox distributions in this region should give better results. 


7. Conclusions 


The numerical method suggested in Section 5 has been applied to 
three examples, which include between them the various edge singularities that 
might be encountered in practice. It has been seen that with a fairly wide mesh 
distribution errors are appreciable only for the edge areas of the mesh. There are 
two ways of avoiding these errors. The first is not to try to calculate the incidence 
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at the singular edge areas, but to extrapolate the values calculated at the inboard 
areas, bearing in mind the types of singularities expected. The second method is to 
decrease the mesh areas locally in the edge regions and so approach the conditions 
at these edges more accurately. However, it is debatable whether the extra work in 
the second method is justified, since it must be realised that the derivation is based 
on a linearised approximation to real conditions, together with the fact that the main 
camber and twist characteristics are given for the wing as a whole to conform to a 
specified loading. 


Appendix 
The numerical method of Example (i) of Section 6 is now discussed in more detail. 


The velocity potential distribution is taken to be 


Where 8=0-6 is assumed. 
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The distribution of 89; ; (0-933)/(Vad) in the mesh pattern chosen for this example is 


shown in Fig. 14, where the incidence is to be calculated for areas inboard from the 
leading edge from Fig. 14(a) and for areas near the leading edge from Fig. 14(b). 


The table of values for A ,;, ;) x, 1 is Shown in Fig. 15. It is presented in the form of a 


V-matrix to ease the computation of 9; ; Ai, 
i, 


The operation of equation (9) for the calculation of y, , can now be applied directly 


and has been applied to the points indicated in Fig. 14. 


From these the average incidence of each of the indicated areas is obtained by 


equations (10) and the results are presented in Fig. 11. 
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Symmetrically Loaded Circular Plates under 
the Combined Action of Lateral and 
End Loading 


RAYMOND HICKS 
(G.E.C. Simon-Carves Atomic Energy Group) 


SumMMaRY: Expressions are obtained for the radial and tangential bending 
moments in a circular plate under the combined action of (a) a lateral load 
concentrated on the circumference of a circle and an end tension or compres- 
sion, and (5) a uniformly distributed lateral load, having a diameter less than 
the diameter of the plate, and an end tension or compression. For both types 
of loading, solutions are obtained for plates which are simply-supported and 
for plates with an arbitrary end rotation. 


In addition, the following limiting cases are considered: (i) concentrated 
lateral load with end tension or compression, and (ii) an infinite plate under 
the simultaneous action of an end tension and a lateral load concentrated on 

the circumference of a circle of finite diameter. 


1. Introduction 


For a laterally loaded plate having small deflections, the bending stress 
distribution is linearly related to the applied load system. That is, any increase in 
the intensity of the applied loading results in a corresponding increase in the 
intensity of stress, with no resulting change in the stress distribution. This implies 
that for a plate under the simultaneous action of two or more lateral load systems, 
the resultant stress distribution can be obtained by adding the stress distribution 
associated with each type of lateral loading. 


For a plate under the combined action of lateral and end loading, the principle 
of superposition is no longer applicable, so that the resultant stress distribution is 
not obtained by simple addition of the stress distributions associated with the lateral 
and end load systems when acting separately. In general the bending stresses due 
to the lateral load are decreased by the application of a tensile end load, and 
increased by the application of a compressive end load. Furthermore, the change in 
the bending stress distribution is not linearly related to the change in the intensity of 
end loading, the relative effect of the end loading being more severe as its intensity 
increases. 
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SYMMETRICALLY LOADED CIRCULAR PLATES 


The deflection and vibration of a circular plate under the simultaneous action 
of lateral and end loading has been considered by Bickley. He confined his 
attention to a plate with a clamped boundary and considered several important 
applied symmetrical load systems as well as a concentrated load at an arbitrary 
point. Bickley was particularly interested in the behaviour of a diaphragm used in 
the construction of a condenser microphone; possibly for this reason he does not 
reduce his solutions to forms which give explicitly the stress distributions when the 
end loading is compressive. 


In the absence of lateral loadings, Bryan’ and Nadai® considered the effect of 


| end compression and determined the load necessary to buckle the plate. They 
| found that the critical compressive stress for a plate with a clamped edge is about 


50 per cent greater than the corresponding stress necessary to buckle a strut having 
a length equal to the diameter of the plate. 


For an end-loaded plate with a hole at its centre, the compressive load necessary 
to cause buckling has been determined by Meissner“. For a plate with a clamped 
outer edge, he found that the buckling stress is a minimum when the diameter of the 
plate is approximately five times the diameter of the hole. 


This paper considers a further series of problems dealing with a symmetrically 
loaded plate under the combined action of lateral and end loading. General 
expressions are obtained for the radial and tangential bending moments and for the 
deflection surface of the plate. It is found that these expressions can be written in 
terms of the modified Bessel functions /,, and K, when the end load is tensile, and in 
terms of the ordinary Bessel functions J, and Y, when the end load is compressive. 
In all cases the order of the Bessel functions is zero or unity, so that numerical 
values for these functions can be obtaiied readily from standard tables, for example, 
Refs. 5, 6, 7 and 8. 


For a plate under the action of a load concentrated on the circumference of a 
circle of arbitrary diameter, the general solution has been used to obtain the moment 
distribution in a circular plate with a concentrated load at its centre. Furthermore, 
this solution has been used to obtain the transcendental equation which defines the 
stress necessary to buckle the plate. Also, the general solution has been used to 
consider the case of an infinite plate under the combined action of an end tension 
and a lateral load concentrated on the circumference of a circle. A solution is also 
obtained for an end-loaded circular plate with a lateral load uniformly distributed 
over a circle of arbitrary diameter; this solution is used to consider the case where 
the load completely covers the plate. For a particular example considered, it is 
found that the effect of an end compression of 4-47 tons/in.? is to increase the 
maximum bending moment by 36 per cent. A further example is used to show the 
attenuation of the bending moment in an infinite plate under the simultaneous 
action of an end tension and a uniform load distributed over a circle of finite 
diameter. 


Finally, where applicable, it is shown that for limiting cases the solutions 
obtained in this paper reduce to those given in Refs. 1, 2 and 3. 
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NOTATION 


w deflection 
r radial co-ordinate 
a radius of lateral load 
b radius of plate 
h_ thickness of plate 
Q, radial shear per unit length 
Q total external shear 
M, radial moment 
M, tangential moment 
L lateral load 
q intensity of lateral loading 
pend load stress 
D flexural rigidity 
vy Poisson’s ratio 
E Young’s modulus 
6 end rotation 
Jn, Yn Bessel functions 
I,,K, modified Bessel functions 
H,® Hankel function 
a= /(ph/D) 
x=rV/a 
~=dw/dr 
C, constant 


2. Fundamental Equation for Small Deflections 


For a symmetrically loaded circular plate having small deflections, the radial 
shear per unit circumferential length is given by the following well known 
expression : — 


d(d’w_ 1 dw 


Referring to the outer section of the plate shown in Figs. 1 and 2, around its inner 
boundary the total external shear acting in a vertical direction is given by 


Q=L-—2taph 
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SYMMETRICALLY LOADED CIRCULAR PLATES 


The negative sign has been inserted to be consistent with the convention that a 
downward deflection is positive and corresponds to a negative radial slope. 


Also, at any arbitrary point in the plate, the external shear is in equilibrium 
with the radial shear Q,. In this case, it follows that 


Combining the expressions for Q,, the equation defining the deflection surface of 
the plate is given by 


a(@w dw) phdw_ 
dr\dr? r ar D dr 2nxrD° 


For convenience let a= Dh 
D 
ph 
dw 
y= 


Using these parameters, the fundamental equation becomes 


Except for the term containing the external shear due to the lateral loading, 
equation (1) is identical to the equation used by Meissner“. 


3. Load Concentrated on the Circumference of a Circle 
(i.e. Concentric Loading) 


For the outer section of the plate shown in Fig. 1, it can be verified that the 
necessary solution of (1) is given by 


_ dw _ 
y= =C,I, (rv 2)+C.K, (rv 2) ‘ (2) 
B 
so that w= a)— Ke a)— log (ri a)— Va (3) 


where B=C,I, (bv 2)—C.K, (bv 2)- 5 log 
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TOTAL LOAD=L 


FicureE 1. Lateral load concentrated on the circumference of a circle. 


Tables and formulae for this Bessel function can be obtained from Refs. 5, 6,7 
and 8. 


Thus, from (2) and (4), 
Cl, (a) +C,K, (ava) - « 


This is the first boundary condition for the outer portion of the plate. 


4. Simply-Supported Plate with Concentric Loading and End Tension 


The radial bending moment is given by the well-known relationship 


dw »vdw 


Substituting (2) and (4) in this expression, two further boundary conditions are 
given by 


hove] + 


, v 
+C, [K, (by a)+ bVa K, (bya) | Ixb?2a3/2D 


270 The Aeronautical Quarterly 


wh 


Us: 


alsc 


whi 


Augu 


(M 
pha— | po 
M M 
a 
b 
cor 
C,: 
wh 
CG 
= 


(4) 


sion 


are 


uarterly 


SYMMETRICALLY LOADED CIRCULAR PLATES 


(M.,=—-DVa {c, [1 (aV2)+ + 


where primes (’) denote differentiation with respect to r. Similarly, for the inner 
portion of the plate, 


dw 


so that [1 . . 
Using the identities 
(z)+ 1, (2) =Zl, (2) (10) 
2K,’ (z)+K, (z)= -2K, (2). 
consideration of boundary conditions gives 
-C 
C= (a/ 2) K, (av 2) K, 2)+ (av 2) SK, (av 
(by 2)+ K, (bv) |} 
1-v) 
C,=-C (bY a)+ dK, (bV2) | I, (a/a)} 
where 


S=K, (av2) [ 1, | + 


also Q=L—2xaphC,I, (av 2). 


Finally, the radial and tangential bending moments are given by equations (12) 
Which follow. 


August 1959 271 


| 
6,7 

(5) 


RAYMOND HICKS 


a)+ 


+ p- (ava) Cul, (av 2) (a<r<b) 


(12) 


Kev 9)|-L (ava) Cul, (av 2) |S (a<r<b) 


Expressions (11) and (12) show that, for a fixed value of the direct stress p, the 
bending moments in the plate are proportional to the lateral concentric loading. 


5. Simply-Supported Plate with a Concentrated Load and End 
Tension 


When the diameter of the concentric loading approaches zero, the solution (12) 
is applicable to a circular plate under the action of a concentrated load L at 
its centre. 


Thus, since =O] (13 


[(avV/ 2) K,(a/2)=I, 


it can be shown that 


ova] - 


r 


Mo= 


1 —v) 
where = (bya (bva) |. 
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K, (2)= -K,’ (2) 


(0-5z)"+?" 


L@= 2 meer) 
y=0-5772. 


Thus, for sufficiently small values of z, 
K, (z)= —log (z/2) 
Using these expressions for K,,(z) and the series for /,,(z), in the limit when p is 
zero, it is found that 
L b 
M,= 4x (1 + ¥) log 


These are well-known expressions for the moments in a laterally loaded plate with 
no end load. 


6. Simply-Supported Plate with a Concentrated Load and End Com- 
pression 


For a plate under the action of a compressive end load, p is negative, so that 
the modified Bessel functions have imaginary arguments. This implies that for a 
compressive end load the radial and tangential bending moments are functions of 
the ordinary Bessel functions J, (z) and Y, (2). 


Now (ze"")=K, (iz)= 7H, 
and A,© (d=J, (z) iT. (z). 
Thus K, (id=- [iJ,(z)+Y, 
Also, it can be shown that 


1, G2z)=J, (2) 


1, (iz)=iJ, (2). 
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Using these relationships, for a compressive end load, it is found that 


+ oy, vay] } (by/2)— ova) | 


It is worth noting that the moments given by equations (15) become infinite when 
bv (2) J, (bv 2)—(1—v) J, (by 2)=0. 


This is the well-known transcendental equation which defines the critical 
compressive stress for an end-loaded circular plate with a simply-supported 
boundary (see, for example, Ref. 9, p. 369). 


7. Arbitrary End Rotation, Concentric Load and End Tension 


This problem is of practical interest, since it has an application in the analysis 
of a cylindrical pressure vessel with a flat head. Also, as a particular example, it 
gives the solution for a plate with clamped edges. 


Referring to equation (3), if # is the radial slope of the plate at r=b, the outer 
boundary condition is given by 


Substituting the appropriate expression for Q, the condition (16) becomes 


a L 
(by a)+C.K, (b/ (av 2) 0. ( 


Proceeding as before it is found that the following expressions are obtained for C,, 
C, and C,. 
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b 
(ava) | G- +7 


Ov aK, (av 2) 2) 41, (av 2) K, (av 2} 


(18) 
(2) 1, (av 2) (b/ (av/a)} 


where 
G={I, (av 2)+ 2) I, (av 2)[K, (av 2) I, (bv 2)—K, (bV 2) 1, 


Using these values for the constants, the radial and tangential bending moments are 
given by (12). For a concentrated load at the centre of the plate, 


—K,(bV2)+ [ (r/a)+ (r/ a) |- 


1, v2) |- 


When p and @ are both zero, these can be reduced to the following well-known 
expressions. | 


For an end compression, these expressions for the bending moments can be 
rewritten in terms of the Bessel functions J, and Y, by using the relationships 
previously given. 


8. Infinite Plate with a Finite Concentric Load and End Tension 


When the stresses acting at the middle plane of the plate are sufficiently large, 
the radial and tangential moments are confined mainly to a region of the plate in 
the vicinity of the concentric lateral loading. This implies that the solution for an 
infinite plate is applicable to a plate of finite size, provided that its least lateral 
dimension is several times the diameter of the concentric loading. It is clear that 
for any particular problem the degree of approximation introduced into the analysis, 
by assuming the plate to be infinite, can be determined readily by comparing the 
stresses and displacements applied to the boundary of the finite plate with the 
stresses and displacements around the corresponding contour in the infinite plate. 
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TOTAL LOAD=L 


FiGurE 2. Lateral load uniformly distributed over a circle. 


Considering (19), if # is zero and b is made infinite, it is found that 


LK, (av z) 


L 


-») 


LK, (av 2) —v) 
L 
(av 2) K, (av 2) [ Tx, 
Pal, 


(20) 


A practical example showing the distribution of the radial moment in an infinite 
plate is considered later. 


9. Uniformly Distributed Load 


Let L be the total load uniformly distributed over the centre portion of the 
plate, as shown in Fig. 2. 


The radial shear at any arbitrary point is given by 
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Thus, for the centre portion of the plate, it can be shown that 


Lr 


v=C,I, (rv 2) 


(22) 


Also, the following expressions give the radial and tangential bending moments 
in the plate: — 


M.=-V@D{ Cl 1 21, 2) | 5K, + 


(av) |G} @<r<b) 


=-V@D{C, 2521, rv 2)- 


(23) 
10. Simply-Supported Plate with a Uniformly Distributed Load and 


End Tension 


For a simply-supported plate, the boundary conditions correspond to the 
following set of equations : — 


Cl, (av 0)+C,K, =0 


(bV/a)+ (bv2) | + Ky (bVa)+ (bv 2) | - 


L(i-y) 
~ xb?a3/2D 


aC, I, (a =0 (24) 


i, (ava) | +C (av 2)+ K, (ay/2) | - 


G-yL 


=0. 


Solving these equations, it is found that C,, C, and C, are given by the following 
expressions. 
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~(aV 2) SK, (aya)+2[s- K, (ava) } 


=c {[x. (by a)+ Yk, |1, (aya) 


(25) 


+ K, (av 2) (21, (av 2)1, (av } 
(1-v) 

(av 12h, (av 2)—(ay2)1, (avy) , 
where C and S are constants defined by (11). 

When a=b, so that the load completely covers the plate, 
- 
- 1 


1, ava] f, 


(1—y) 
where 1, (ava)— 1, (aV/a)|. 


Substituting the appropriate series for the Bessel functions, it can be shown that 
when p is zero the bending moments are given by the following well-known 
expressions : — 


M,=q 
q 3 
where q is the intensity of the lateral loading. 


For many practical problems the parameter (a¥/ 2) is appreciably less than two, 
in this case it is convenient to write the expressions (26) in a series form. 
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- 
Now /7/,(z)=1+ (mle 


27"+1 (m! 1)! 


1, @=5+. > 


so that 
alll tpl msi 
Substituting these series in (26), it is found that 


{ 


(3+) (1+ 39)] 


22"+1 


My= {14 
x 


1 
where (m!)? (m+ 1) 


and M,, and Mg, are the corresponding moments in a plate with no end load. 


Provided that (a\/ 2) is reasonably small, the bending moments can be evaluated 
by considering the first few terms in the series. For example, suppose the end 
tension is 4-47 tons/in.?, E= 13,400 tons/in.*, v=0-3, and a=15h. Then (a/2)~1 
and the maximum radial bending moment is given with sufficient accuracy by taking 
the first term in the series. Thus if Minax. is the maximum moment in a similarly 
loaded plate with no end load, then 


1+0:067a?2+. . 

where the end load is tensile; this gives (M,) max =0°81M max o- 
Similarly, when the end load is compressive, (M,)max= 1°36M max o- 


It can be verified that, for a plate with end tension, the maximum moment acts in a 
tadial direction; while, for a plate with end compression, the maximum moment acts 
in a tangential direction. For the particular example considered, the maximum 
tangential moment with end compression is given by 


(Mo) max = 1-43 
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11. Arbitrary End Rotation, Uniform Distributed Load, and Eng 


Tension 


For an arbitrary rotation of the plate at r=b, the boundary conditions 


correspond to the following set of equations : — 


(aV CK, (ay =0 


Cy 2)+ + 2 (a2) — 9-0 

C, a)+ ave” 2) | +C, (av 2)+ = 
L 


(ava)+ =(), 


= 
Solving these equations, 


bé 
+ 


| 


L 
2GrazDI, (av 2) 


L 
~ (av 2) 


(bv 2) K, (av 2) 1, (av 2)—K, (av 2) [21, (aV'2)- 


a) 1, 2)]+ 2 61, (av 2) K, } 


LI, (av 2) 


2GrazD 


[2/, (av z)—(av 2) (av 2)], 
where G is given by (18). 


For a load uniformly distributed over the whole area of the plate, 


_ 


280 


(27) 


(28) 


+ (29) 
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Ficure 3. Example of a radial moment attenuation in an infinite plate. 


When p and @ are zero, it can be verified that 


where q is the intensity of loading. 


12. Infinite Plate with End Tension and a Uniformly Distributed 


Circular Load 


For a load distributed over a finite circular area having a diameter 2a, the 


limiting values of (28) can be reduced to 


CHL [2K (av z)—(aV 2)K, (av 2)] 
(av 2) K, (av 2) 
C,=0 


L 
2nazDK, (av 2) 


Substituting these values for the constants in (23), the resulting expressions give 


the radial and tangential bending moments in the infinite 


plate. 


When the diameter of the applied loading is small and the stress at the middle 
plane of the plate is reasonably high, the arguments of the Bessel functions become 


small compared with unity. 
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In this case reduction of expressions (23) and (30) gives 
(1+) 
M,=Me=~ -L- log ) (0<r<a) 


ra 


L 


= 


Also, for small values of r > a, it can be shown that 


To illustrate the attenuation of the radial bending moment in an infinite plate, 
consider the example where the diameter of the loading is three times the thickness 
of the plate, p=4-47 tons/in.?, E=13,400 tons/in.?, and v=0-3. Thus, in this 
example a/2=0-1, so that the approximate expressions (31) can be used to 
determine the moment distribution. This is plotted in Fig. 3. 
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The Drag of Ducted Bodies with Annular 
or Side Intakes in Supersonic Flow’ 


H. PORTNOY 
(Royal Technical College, Salford) 


SuMMaARY: A slender ducted body of arbitrary cross section is studied. This 
possesses an annular, side or “scoop” type of intake at some streamwise station. 
The parts of the body forward of and behind the intake are both permitted 
to have discontinuities of longitudinal profile slope at streamwise stations 
widely separated both from each other and from the intake. The external 
flow is taken to be entirely supersonic and it is assumed that the intake lip 
is sharp and that there is no “spill over” due to choking of the internal flow, 
which must be supersonic for a short distance inside the duct. 


The interference drag between the forebody and afterbody is calculated 
and, since their isolated drags are known from previous work, the drag of the 
whole body is then known. 


The results are applied to bodies of elliptic cross section and some 
numerical comparisons are made between different bodies of this type. 


The work of this paper is applicable only to intakes with sharp lips of 
small external streamwise slope (because of the slenderness requirement) under 
conditions of no spillage. 


l. Introduction 


Ducted slender bodies in supersonic flow have been considered by several 
authors. Lighthill® has dealt with bodies of revolution whose meridian section 
may have discontinuities of slope, and a ducted body of revolution with an open 
nose is a special case of these results. In a similar way the results of Fraenkel® 
cover the case of an open-nosed ducted body of elliptic cross section, the eccentricity 
being constant, while more general open-nosed elliptic bodies have been dealt with 
in Ref. 3. Finally, general open-nosed slender bodies of any cross section form a 
special case of the theory of Ref. 4. In Ref. 4 it was indicated how the results could 
be applied to other types of ducted body, such as bodies with annular or side intakes, 
but they were not explicitly applied to the drag problem for such configurations, 


*This paper is extracted from a thesis for a Ph.D. degree in the University of London, written 
while the author was at the Imperial College of Science and Technology. 
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although the external side force and moment expressions were derived. It is the 
purpose of this present paper, using a result of Ref. 4, to deal with the external drag 
problem. For the particular case of a ducted body of revolution with an annular 
intake, a complete set of results for the external forces, together with an expression 
for the internal lift force, was derived by Ward®. In all the other references only 
the forces due to external pressures were considered. 


The type of configuration considered here consists of a slender forebody, in 
which for generality we permit a finite number of widely spaced discontinuities of 
longitudinal profile slope, after the manner of Ref. 4, followed by a ducted afterbody, 
which, similarly, may have discontinuities which are far from the intake and are 
widely separated from each other. This ducted afterbody may have an intake of 
annular type (not necessarily circular) or it may have side intakes or a “scoop” 
intake. Whichever configuration is considered it is assumed that the external flow 
is independent of the internal flow; which implies that there is no “spill-over” due 
to choking of the internal flow, which must be supersonic for a short distance inside 
the body. 


As an exaniple of the method, some results and comparisons for the case of a 
body of elliptic cross section are given. 


NOTATION 

A(z) major radius of an elliptic body cross section at station z 

b, see equation (4) 
boo, box +=see equation (6) 

B(z) minor radius of an elliptic body cross section at station z 
C contour of a section z=constant on body 
D drag. Various suffixes define the component in question 
D, interference drag between forebody and afterbody 


e; values of z at which the longitudinal surface slope of a slender body 
is discontinuous 


F, (S’, 8) see equation (8) 
F,(S’,8) see equation (10) 
H(z) Heaviside unit function. H=1 when z>0; H=0 when z<0 


k, a quantity defined as 1 when considering the drag of an afterbody 
in the field of a forebody and as 0 when the afterbody is isolated 


1 length of a body 
m_ number of discontinuities in an afterbody 
M free stream Mach number 
r,9,z cylindrical polar co-ordinates 
S(z) cross-sectional area of a slender body at station z 


t the complex variable Ae‘# 
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V, free stream velocity 
x,y,z Cartesian co-ordinates 
Z, see equation (13) 
B= /(M?-1) 
See equation (18) 
A symbol indicating a jump or change in a quantity 
€ slender-body thickness ratio. A slenderness parameter 
¢ the complex variable x+iy 
A» =(A + B)/2 
y,7 normal and circumferential variables on a contour z=constant 
p, free stream density 
$%, two-dimensional harmonic potential 
w see equation (11) 


Suffixes 
i,i+1 relates to a segment of slender body between the i and the (i+ 1)" 
discontinuity 


a,b,c _ relate to the bodies a, b and c. 


2. The Drag of a Segment of a Slender Body Between 
Two Discontinuities 


Let Ox, Oy and Oz be a right-handed system of axes with Oz lying in the 
direction of the undisturbed stream and the origin at the nose of the forebody*. 


A typical body segment is considered, lying between the streamwise stations 
z=e; and z=e;,,, at which not only the body streamwise slope but also the body 
profile radius may be discontinuous. However, as high speed exhaust gases cannot 
in general be treated by the method of small perturbations which is used here, a 
segment influenced by such gases cannot be dealt with, and so in practice the 
discontinuities of body radius must be such that the radius increases with increasing 
z, Le. we permit only intake openings, except at the base, which could coincide with 
Z=é,,. However, it is convenient in writing expressions for wave drag to assume 
the body continued downstream of the base by an infinite cylindrical extension with 
the base cross section. This extension does not influence the wave drag upstream of 
the base in supersonic flow, but enables one to define slope discontinuities at the 
base and so to write wave-drag expressions in a concise manner. The normal to a 
body contour in a plane z=constant is denoted by v and the circumferential variable 
is represented by 7, while the body cross-sectional area at station z is written S (2). 
The slenderness parameter ¢ may be taken as the greatest value of the streamwise 
slope of the body and, for the body to be slender, S (z)=O (¢*) everywhere. 


*In the most general case admissible under the present analysis, i.e. that where the forebody 
has an open nose, the origin is taken at the centroid of the nose section. 
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In Ref. 4 it is shown that the drag, D,, :,,, of the typical segment is given by 


| 8” (2) b, @ dz+S' by —) b, -) | +0). 


where ¢, is a two-dimensional harmonic potential in a plane z=constant which 
satisfies the conditions 


dv 
dz 
on the body contour, and 
(z) 
lim “ie log r) =0, (3) 
r being equal to /(x?+ y’) and 
by tog log (z—z,) dS’ (z,). 


C (i4:)- denotes the body contour just upstream of z=e;,, and so on. In equation 
(4), B= /(M?-1) and the integral is a Stieltjes integral. 


3. The Drag of a Body with Side or Annular Intakes. Interference 
Between the Forebody and Afterbody 


The result (1) is now used to derive the interference drag between the forebody 
and afterbody of a configuration with intakes at a single streamwise station as 
described in the Introduction (see also Fig. 1). The drags of the forebody and 
afterbody being known from Ref. 4, this is the only unknown component of 
the drag. 


The forebody is unaffected by the afterbody in supersonic flow, and so the 
interference drag is the drag of the afterbody in the presence of the forebody minus 
its drag when isolated (at the same free stream conditions). 


The afterbody has discontinuities at z=e;(i=0, 1,2. . . m), where z=e, is the 
position of the intake and z=e,, is the position of the base. 


Thus, if base drag is ignored, the following expression for the afterbody wav 
drag may be obtained by summation of (1). This expression covers both the cast 
where the forebody is present and when it is absent. 
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by 


FOREBODY~ 


m “DISCONTINUITIES” IN AFTERBODY 


which FiGurE 1. 
() Wave drag of afterbody § [ dv 
4p,.V,? ? 2k; dr+ 
egt 
dv dy 


+2 [ | Ss” (z) b, (z) dz+ b, (e:—) AS’; + S’ (e,+) b, (e, -)| +O (>), (5) 
(4 ‘ 


where k; is defined as unity when the forebody effect is to be included and as zero 


lation | otherwise. b, must, of course, also be defined according to the conditions. The 
usual notation 
Afi=f (x, y,e:+)—f (x, y.e:—) 
ence 
isused. Let 
*body 
On as 1 B 
bu= {s@ log H (z—e0)— | log (z—z,) dS’ (z,)—S’ (e, +) log (z—e,) 
nt of ‘ “oF 6 
1 
atte | 8’ log H (e,—2)- | log (@—z,) dS" 
ninus 
. where H is the Heaviside unit function and log (z—e,) is taken as zero for z < @,. 
s tht | It is evident that, for the afterbody 
wave Des, when the forebody is absent; 
Case 


1, when the forebody is present. 
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Hence, if D; is the required interference drag, subtraction yields 


_Dy =4 dv 
+2 [ | Ss” (z) box (z) dz+ = Dox (e,—) AS’, + S’ (e,+ ) Dox (e, -) +0 
= (2) dr+F, (S’, B)+O(€ 
\Gz Janes . . . . (1) 
where 


F,(’, 8)= | S” (2) [s (e, —) log (z—e,) | log (z—z,) dS’ dz— 


£9 o-— 


m 


[as’ tog (e,—z,) dS" (z,)- AS’, (e, -) log (e -e,)| + 


+S’ (e,+) (e, —) log ~ | log (e, —z,) dS’ (2, | 
== J toe 377,48 (z,) dS’ (z)+ S’ (e, +) | log (z)+ 


ent 
+S’ (e,-)log —S'(e,-) | as’ (2 | @ 


The drag result of Ref. 5 is readily recovered from equation (7), making use of 
expressions for the drag of the forebody and afterbody obtained from Ref. | (for 
example) and the fact that, for a body of revolution, 


_ 
log r. 


Equations (7) and (8) also form an extension of equation (8.14) of Ref. 2, which 
gives the interference drag of an elliptic cross-sectioned configuration in which the 
forebody and afterbody are separated by a plane ‘~ which the longitudinal profile 
slope only is discontinuous (and not the body radiv »). 
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4, The Afterbody of Elliptic Cross Section 


Next consider the case of a smooth afterbody having two longitudinal planes of 
symmetry and elliptic cross sections. For the present, no particular restrictions, 
apart from slenderness, will be imposed on the forebody. 


The limiting case in which the elliptic intake opening is tangential to the fore- 
body at two points is of some practical interest, as then there is something 
approximating to a body with side openings. 


The drag of the isolated afterbody may be written down directly from the 
general result derived in Ref. 3, equations (10) and (11). The only “discontinuities” 
of the afterbody, in this case, occur at the intake end (z=e,) and the base (z=e,); 
therefore 


e€ e 
1 


| toz _ (z,) dzdz, + 5 tos (z) dz— 
toV, 2x | z-2 | 


(e,-) log S” (z) dz—-S’ (e, +) (e, —) log \ 


where A (z) and B(z) are the major and minor radii of the section at station z, and 


1 
F, (S’, 8)= | Iz (z) S” (Z, ) dzdz,+ — (e,+) {log - (2) dz— 
£9 
—_ 
-S'(¢,-) | log (2) 
"0 


+5 {is @, } log (10) 
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To find the interference drag D, it is first necessary to evaluate the contour 
integral appearing in (7). The body section at z on the after portion will, in this 
case, be given by the parametric equations 


For a body given in such a parametric form it is readily shown that, if ) 
increases anticlockwise, 


dy = 


In the present instance 


= [5 +48 cos 24] sinty+ 


Making the conformal transformation 


(13) 


x+iy=(=t+ [A? (z,)—B? (z,)]/(40) 


where t=Ae/#, the body section at z=z becomes a circle of radius A=A,=(A + B)/2. 
Also on the body » =~ and 


The boundary condition (2) for the potential of the afterbody may now be written as 


_ dt _ dy 
Ov body ~ dz 
On Ss’ A 
and so = +4B? (z cos 2u} 
on the afterbody. 
Now since 
the contour integral becomes 
20 3 
0 
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But, as the forebody contour lies within the ellipse at z=e,+, the transformed 
forebody contour will lie entirely within the circle A=A,, as do the singularities of 
the transformation, and hence so will all the singularities of (%,):,--,-. Therefore 
this potential may be expanded on the circle A=A, as a series in the form 


2 =(44B)/2 


m=1 


where Y¥m=Ymr+iym1, and the coefficient of log ¢ has been immediately placed equal 
to S’(e,—)/2%, as this is readily derivable using the expression analogous to (16) 
for the contour C-,-, together with Cauchy’s theorem. 


On inserting (15) and (18) in (17) one obtains 


\dz 


Equation (7) now gives, upon substituting this result, 


(0) 


It is worth noting two results which follow from (20). 


(i) If the eccentricity of the afterbody sections is constant, i.e. d (A/B)/dz=0, then 
D,; does not depend on the cross-sectional shape of the forebody, but only 
on S’ (e, —). 


(ii) If the forebody is circular, 
S’ (€.—) 
=F . 


+ F, (S’, 8)+ O 
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5. Example 


The foregoing results are now used to compare the drags of three ducted bodies 
having the same distribution of cross-sectional area, the same forebody, and 
possessing two streamwise planes of symmetry. The forebody is a circular cone of 
semi-angle «, and it extends to half way along the total body length, this last 
quantity being denoted by /. The ratio S (e,+)/S(e, —) is taken to be 1°5. 


The three bodies differ in the shape of the after portion. The first (body a) 
has an axi-symmetric afterbody in the form of a cone frustum also of semi-angle <, 
The second configuration (body 5b) has elliptic afterbody sections varying from a 
section tangential to the forebody at top and bottom, at the intake, to a circular 
base, d(A/B)/dz being held constant. Finally, body c also has elliptic sections and 
an intake opening tangential to the forebody at the top and bottom, but its sections 
have constant eccentricity. Their drags are respectively denoted by D,, D,, and D.. 
It is found in the Appendix, from equations (9) and (21), that 


and 

the errors being O (:°), and, using the result of Ref. 5 to obtain D, (see Appendix), 


These last two expressions are plotted in Fig. 2 against 2:. It should be noted that 
the drag differences shown in (22) and (23) are independent of Mach number, as 
this only enters the drag expressions through the functions F, (S’, 8) and F, (S’, 8) of 
equations (21) and (9), respectively, which are eliminated in subtraction, as all the 
bodies have the same S’ (z) distribution. 


Evidently body c has the smallest drag, followed by a and lastly by b. This is 
similar to the result of Ref. 2, where it is shown that open-nosed elliptic bodies of 
constant eccentricity have less drag than the body of revolution with the same cross- 
sectional area distribution. It appears that changing to elliptic sections tends to 
decrease the drag, but varying the eccentricities of the sections has an adverse effect 
which can completely overcome the advantage of the elliptic sections, as here. 


From Fig. 2 it is seen that neither b nor c has a drag more than a few per cent 
different from that of the equivalent body of revolution a, in the range considered 
(which is a reasonable working range for the theory). Thus it may be expected that 
the tangential elliptic opening is as good an approximation to any of the usual side 
opening intakes as will ever be required in practice, since these are generally closet 
to the tangential elliptic opening than it is itself to the circular intake. 
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0:01 
0 0-02 003 006 008 O10 O@ OF O16 0-20 
ALL THREE BODIES HAVE THE SAME 
-001 (z) DISTRIBUTION AND THE SAME FOREBODY __| 
gopya 
+ BODY OF REVOLUTION 
‘A 
-0-02 (2) =cowstant 
AFTERBODY ELLIPTIC 
AT INTAKE CHANGING 
TO CIRCULAR BASE 
io WITH CONSTANT 
ECCENTRICITY 
-0'04 Og 
-0:05 
-0°06 


FIGURE 2. 


As D, and D, differ so little from D,, Ward’s result®? would be adequate 
practically for finding these drags except for very low values of ¢, especially since 
the differences are of the same order, or approaching the same order, as the error 
terms, except for these low values. 


However, the configurations b and c are not very different from bodies of 
revolution and it is to be expected that for bodies diverging further from the axi- 
symmetric shape, e.g. much flatter elliptic sections with larger variations of 
eccentricity, or bodies carrying nacelles on narrow wings (which could be treated by 
the present general method if certain restrictions were fulfilled), the differences from 
the axi-symmetric result would not be negligible. This would also, most probably, 
be the case for bodies not possessing two planes of symmetry, such as bodies with 
“scoop” intakes. 
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Appendix 
All three bodies described in the example in Section 5 have the same area distribu- 
tion, given by 


S(z)=7¢?2", 0<z</1/2 } (26) 
and { 1/2(/ 1-5-1) +z } 2, l/2<z<l 
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The station z=e, is therefore at z=//2. 


In addition, for bodies b and c, the value of the ratio A/B for the elliptic intake 
opening is seen to be 1-5. For body c this value of A/B is maintained until the base, 
whereas for body 5 it varies at a constant rate to produce a circular base. Thus the 
constant value of d(A /B)/dz equals —1/1. 


Since all three bodies have the same axi-symmetric forebody, the difference of 
the drags of any two of them equals the sum of the differences of their isolated 
afterbody drags, as given by equation (9), and their D,’s, as given by (21). Also, 
since S’(z) is the same for all three, the functions F, and F, disappear on subtraction as 
they depend only on S’(z) and ~. 


Thus, using (26) with equations (9) and (21), 


D,—D, (og el 


(we? 1) (4/ 1-5 me? I) ViSel_ 1-25 el 


($3) (F)+0@ 


which is equation (22). 


+p 


In the same way equation (23) is derived. 


The result of Ref. 5 is for the drag of the whole of a ducted body of revolution of 
unit length without discontinuities other than those at the intake station and the base. 
In the present notation, and adjusted to a length /, it reads 


{tog Se de, ae FO 5”(2) d+ 
0 0 


S” (z) dz— 


[S’()}* S’ (e,+)—S’ (e,—) 


Fle,+) 


leg—zl 


+0 (25) . (28) 
The error term is that found in Ref. 4, as used throughout this paper, whereas Ward gives 
O (e® log? «). 
Equation (28) can be obtained from the results of the present paper as a special case. 
Substitution of (26) in (28) yields 
D, 
tp, V. 
This result, together with (22) and (23), gives (24) and (25). 


e4 (2:50 log +0:548)+O(e5). . & 
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Review 


The Potential Theory of Unsteady Supersonic Flow. J. W. Miles. Cambridge 
University Press, London, 1959. 220 pp. Illustrated. 45s. 


“Observe how ... the critic is invariably a tremendous specialist in the subject under 
review, and must at all costs be more so than the author of the book discussed. It 
doesn’t matter if the subject is as remote as the study of Greek in Lower California, the 
reviewer must be there before.’’’) It seems only fair to warn readers at the outset that 
this fundamental principle will not be upheld in the present review, for the reviewer's 
knowledge of unsteady supersonic flow, even though it were fortified by every gambit in 
the Lifeman’s armoury, would still be no match for Professor Miles’s complete mastery 
of his subject. 


In a sense the problems of unsteady linearised wing theory are of a classical type, 
and may be viewed as an exercise (for advanced students) in the theory of sound: if the 
co-ordinate system is fixed in still air, one seeks solutions of the wave equation which 
satisfy certain mixed boundary conditions in the plane of the wing, and this may be 
achieved, in principle, by a suitable distribution of retarded-potential doublets in that 
plane. Unfortunately this approach (due to Kirchhoff'?’) leads, in general, to rather 
terrible integral equations for the doublet strength (although for points influenced only 
by a supersonic leading edge the problem is simpler: sources of known strength may be 
used), and modern workers have tended to explore a number of other avenues. It is 
this recent work which Professor Miles expounds, and one is left with the following 
general impression. 


The choice of a best co-ordinate system is not obvious: co-ordinates fixed in still 
air are “natural” for the potential equation and the pressure relation, and co-ordinates 
moving with the wing for the boundary condition: other systems, with unnatural times, 
may effect a compromise. Integral transforms (one, two, or three may be used, 
according to the problem) are invaluable in reducing the algebra to reasonable propor- 
tions. The transformed problem is straightforward for two-dimensional wings and 
those with supersonic leading edges; it turns out to be an old one in diffraction theory 
for the semi-infinite rectangular wing; and it has a solution in Mathieu functions for the 
complete rectangular wing. The inversion and numerical evaluation of this last solution 
is not without difficulties, and for wings with subsonic side edges an extension of Evvard’s 
method offers the best hope of success. But for the delta wing with subsonic edges one 
must accept defeat, and expand in powers of the frequency. Finally, if the aspect ratio 
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can possibly be considered “small”, one is in the calmer realm of slender-wing theory, 
but even here high frequencies may require the addition of a time derivative to the 
transverse-plane Laplace equation, and thus make the problem more difficult than in 
steady flow. 


Professor Miles steps nimbly through this material, with as much unity of method 
as the vagaries of the subject allow. In addition to the acoustic and slender-wing 
theories already mentioned, he systematically derives six other approximate equations 
for the velocity potential, together with the restrictions on reduced frequency, amplitude, 
aspect ratio, and Mach number appropriate to each. Slender non-planar bodies are 
also treated, non-linear problems are described briefly, and an appendix deals with reverse 
flow theorems. The main applications are to damping in pitch and gust loading. There 
is no space for such questions as accelerated flight, drag, or the calculation of all the 
relevant forces and derivatives, but ample references are given to papers on those topics. 


The reviewer was able to detect only two minor flaws. First, the restriction toa 
perfect gas (p. 1) is unnecessary for the linearised approximations with which the book 
is mainly concerned. Secondly, the special notation dn/0x and dn/dt of slender-body 
theory (pp. 158-9) should be defined, for the function n(x,s,t), which the author uses 
explicitly, is meaningless and does not exist in general. 


The reader is expected to concentrate and to skip no pages; constant reference is 
made to the material in earlier chapters (usually by equation numbers only), and the 
pace of the mathematics may leave the reader somewhat breathless at first—but 
gradually a definite pattern emerges. 


This authoritative account of an important subject, with its complete guide to the 
existing literature, is a most welcome addition to the aerodynamicist’s library. 


L. E. FRAENKEL. 
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